Theory of Cold Atoms: Basics of Quantum Statistics 

V.I. Yukalov 

en ■ Bogoluhov Laboratory of Theoretical Physics, 



O 
(N 

< 

in 

(N 



C/3 



X 



Joint Institute for Nuclear Research, Dubna 141980, Russia 



Abstract 



i-C . The aim of this Tutorial is to present the basic mathematical techniques required for an 

Q ! accurate description of cold trapped atoms, both Bose and Fermi. The term cold implies that 
^ I considered temperatures are low, such that quantum theory is necessary, even if temperatures 
are finite. And the term atoms means that the considered particles are structureless, being 
defined by their masses and mutual interactions. Atoms are trapped in the sense that they 
form a finite quantum system, though their number can be very large allowing for the use 
ci , of the methods of statistical mechanics. This Tutorial is the first part of several tutorials, 
Vh ! giving general mathematical techniques for both types of particle statistics. The following 
'-^ I tutorials will be devoted separately to Bose atoms and Fermi atoms. The necessity of 
C ■ carefully explaining basic techniques is important for avoiding numerous misconceptions 
Q ■ often propagating in literature. 
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1 Introduction 

Physics of cold atoms has recently been attracting high attention because of the breakthrough 
in experimental techniques of cooling and trapping atoms, resulting in a series of experiments 
with cold trapped atoms, when quantum effects become principal. In these experiments, 
different atoms can be prepared and studied in a number of states, with varying density, 
temperature, atomic interactions, and external fields. This allows for a detailed comparison 
of theory and experiment. However, in order that such a comparison would be meaningful, 
the used theory must be self-consistent and mathematically correct. In the other case, 
various misconceptions and paradoxes can arise. In theoretical literature, one can meet 
numerous incorrect statements and unjustified approximations, which especially concerns 
systems with Bose-Einstein condensate. It looks, therefore, that a Tutorial, explaining in 
detail the mathematical techniques that are necessary for a correct description of cold atoms 
is timely and important. 

The present Tutorial is devoted to the general techniques of characterizing atomic systems 
for which quantum theory is necessary. Applications to Bose and Fermi systems are planned 
to be treated separately in the following tutorials. To make the presentation self-consistent, 
it is written with sufficient details, in order that the reader would not be forced to jump 
to other numerous publications. The related literature is cited. But, since the Tutorial is 
not a review, the citations are mainly limited to books, review articles, and most important 
original papers. Some topics, relevant to the physics of cold atoms have been described in 
the books [1^1] and review articles [51-122]. Several results, presented in this Tutorial are new. 

Throughout the text, the system of units is employed, where the Planck constant /i — )■ 1 
and the Boltzmann constant ks -^ ^■ 

2 Phenomenological relations 

Before starting the exposition of a rigorous microscopic approach, it is useful to recall in brief 
the main phenomenological relations that often are common for quantum as well as classical 
systems and, when applied to quantum systems, describe the relations between their average 
characteristics. 

2.1 Thermodynamic quantities 

Thermodynamic relations are the same for either quantum or classical systems, provided 
that the considered system is in equilibrium and is sufficiently large, with its number of 
atoms N and volume V such that A^ ^ 1 and the volume being proportional to A^. 
The transition to bulk matter corresponds to the thermodynamic lim,it 

N 
A^ — )■ oo, V^ — )■ oo, — — )■ const , (1) 

so that the average atomic density p = N/V is finite, < p < oo. 

In the case of finite systems, it is not necessary to invoke the thermodynamic limit, 
though, for the validity of thermodynamic relations, the system is to be sufficiently large, as 
is mentioned above. When atoms are trapped in an external confining potential, instead of 



limit (1), one defines [23] the effective thermodynamic lim,it 

A^ — )■ oo, A]s[ — )■ oo, — — )■ const , (2) 

where A]^ is an extensive observable quantity. If A]^ is proportional to V , then the effective 
limit (2) reduces to the standard limit (1). 

The notion of thermodynamic limit is important not only for bulk matter, but for finite 
systems as well, being closely related to system stability. Even dealing with finite systems, 
it is instructive to check the validity of thermodynamic relations under the effective thermo- 
dynamic limit, which makes it possible to understand whether the finite system is stable. 

Thermodynamic quantities are distinguished into intensive and extensive. Intensive 
quantities, being defined for a large system, do not depend on asymptotically large A^, 
remaining finite for N -^ oo. Examples of such intensive quantities are temperature T, pres- 
sure P, chemical potential /i, magnetic field B, electric field E, specific heat Cy, isothermal 
compressibihty kti and other susceptibilities. 

Contrary to this, extensive quantities are proportional to A^ in some positive power, so 
that they tend to infinity together with A^. For systems confined in a volume V , extensive 
quantities are proportional to V as well as to A^. When the atomic system is trapped 
in a confining potential, and the volume is not defined, then extensive quantities can be 
proportional to A^", with a < 1. Examples of extensive quantities are energy E, entropy S", 
magnetic moment (magnetization) M, electric moment (polarization) P, and thermodynamic 
potentials. 

Thermodynamic potentials define the thermodynamic state of the system. Each thermo- 
dynamic potential is a function of those thermodynamic variables that are necessary and 
sufficient for a unique description of thermodynamics of an equilibrium system. The most 
often used thermodynamic potentials are as follows. 

Internal energy E = E{S, V, N, M). In the presence of an electric field E, the latter enters 
all equations similarly to the magnetic field B, because of which we show below only the 
dependence on B. In the same way, electric polarization P enters similarly to magnetization 
M. The internal energy differential is 



Hence, one has 



dE = TdS - PdV + fxdN + B-dM. (3) 



dS ' dV ^ dN ' 9M" ' ^ ' 



where a = x,y, z. 

Free energy F = F{T, V, N, M). It is related to the internal energy as 

F = E-TS , (5) 

and its differential is 

dF = -SdT - PdV + ^dN + B • rfM . (6) 

Therefore, one has 



ar ' dV "^ dN ' 9M" 



Gibbs potential G = G{T, P, N, B). It can be represented as 

G = F + PV -B-M = fiN . (7) 

And its differential is 

dG = -SdT + VdP + fidN -M dB (8) 

From here it follows: 



dT ' dP '^ dN ' 95" ■ 

Grand potential Vt = VL{T, V, /i, M). It can be written as 

n = F- fiN = -PV + B ■ M , (9) 

with the differential 

dVt = -SdT - PdV - Ndfi + BdM. (10) 

Respectively, 

s=_e. p=-^ iv=-a, B" ^" 



dT ' dV 9/i ' aM" ■ 

Some other thermodynamic potentials can be found in Kubo [23]. 

2.2 Equilibrium and stability 

An equilibrium state of a system corresponds to an extremum of a thermodynamic potential, 
when 

SE = 6F = 6G = 5n = 0, (11) 

with respect to the variation of parameters. 

An equilibrium can be stable, metastable, or unstable. The absolute minimum of a 
thermodynamic potential corresponds to a thermodynamically absolutely stable system. A 
relative minimum represents a metastable system. And a maximum characterizes an unstable 
system. The sufficient conditions for the thermodynamic stability are 

6^E>0, 6^F>0, 6^G>0, 6^Q > . (12) 

These stability conditions can be expressed through second derivatives of thermodynamic 
potentials defining the related susceptibilities. 

The system response to temperature variation is described by specific heat (heat capacity) 



Since S = —{dF/dT)v, we have 

N \dT 



c^-im)^- M 



Using the relations F = 


= E- 


- TS and 

WJv 


we get 




Cv^l 



\dTyv 
The isothermal compressibility describes the system response to pressure variation, 



Substituting V = {dG/dP)TN gives 



1 fd^G\ 



The isothermal compressibility can also be written as 

1 fdP\'^ 1 fd'^F 

Kt = 



v\dv)^^ v\dv\^^ 



And, with p = N/V , one has 



.,.1(1^) ^U^Y . (19, 

The isothermal compressibility is a very important quantity characterizing the system 
thermodynamic stability, because of which it is useful to have its other representations. Let 
there be no magnetic and electric fields. Then the Gibbs potential differential is dG = 
—SdT + VdP + pdN . Since G = fiN, it follows Ndfi = —SdT + VdP, from where 





TV 


1 fdP^ 

" P \dN, 


' TV 


ON' 


'TV 


rdN\ 

^'[dpj 


TV 



Inverting this yields 



For the related variables N,P,V, one can write 

Under the fixed N, implying the zero total variation dN = 0, we have 

dN\ rdP\ fdN\ 

'dp)y\dv)^^[dv)^-^' 



from where 



Using the Maxwell relations 



and the equalities 



we get 






rdV\ fdfi 



ydNj^p \dP 



p 



TN 



Hence 



This gives 






9P J TV V ^-P / TN 



[dpj^j, ~~~p [dpj^y ~'7 V^y 



TV 



''^=N\dp)^y~Wp\'djl)^y-'^\W)Tv' 

The isothermal compressibility is closely connected with the sound velocity s, defined by 
the equation 



where m is atomic mass. As is clear, 



«r = -^ . (22) 



Magnetic susceptibility is defined through the susceptibility tensor 



1 /(9M"\ 

Using here M" = — dG/dB°' gives 

- _ 1 f^!^\ 
^-/3 ^ [dB-dBf^J ■ ^ ^ 

Diagonal parts Xaa of tensor (23) are called magnetic susceptibilities. 

The thermodynamic characteristics, describing the system response to the variation of 
thermodynamic variables and which are expressed through the second derivatives of ther- 
modynamic potentials, all can be termed susceptibilities. Thus, the specific heat Cy is 
thermal susceptibility (susceptibility with respect to temperature variation), the isothermal 
compressibility kt is mechanical susceptibility (susceptibility with respect to pressure vari- 
ation), and the magnetic susceptibility Xaa defines the system response to magnetic field 
variation. According to the stability conditions (12), all these susceptibilities are to be 
non-negative. Being intensive quantities, the susceptibilities of stable equilibrium systems 
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have to be finite, even in thermodynamic hmit. In this way, the necessary conditions for 
thermodynamic stability of an equihbrium system are: 

0<Cy<OO, 0<fi;T<OO, 0< Xaa < OO . (25) 

These stabihty conditions are the direct consequence of the convexity properties (12) of 
thermodynamic potentials and of the fact that susceptibihties are intensive quantities. The 
compulsory finiteness of the susceptibilities is easily understood. Take, for instance, specific 
heat Cv- If it would be infinite, then, according to Eq. (15), an infinitesimally weak tem- 
perature increase would yield infinite system energy, which is meaningless. If compressibility 
kt would be infinite, this, by definition (16), would imply that an infinitesimally small fiuc- 
tuation of pressure would make the system volume either zero or infinity. Similar arguments 
are applicable to magnetic susceptibility. 

2.3 Phases and symmetries 

Thermodynamic phase, defined by a set of thermodynamic variables, is a system state qual- 
itatively different from other states of the same system [23]. To make this vague definition 
more precise, one introduces order parameters such that to distinct phases there would corre- 
spond different orders. An order parameter is a quantity specifying a given thermodynamic 
phase and qualitatively distinguishing it from other phases. 

For example, the distinction between ferromagnetic and paramagnetic phases is done by 
means of magnetization M, so that M = |M| is nonzero for ferromagnetic phase (M > 0, 
ferromagnet) and is zero for paramagnetic phase (M = 0, paramagnet). Similarly, the 
nonzero or zero polarization P distinguishes ferroelectric from paraelectric phases. 

CrystaUine phase is characterized by a set of spatial vectors {a}, called lattice vectors, 
with respect to which the particle density is periodic: p(r -|- a) = p(r). While liquid-phase 
density is uniform: p(r) = p. The role of an order parameter, distinguishing between the 
crystalline and liquid phases can be played by the difference 

Ap = maxp(r) — minp(r) . 

r r 

For crystalline phase, Ap > 0, while for liquid phase, Ap = 0. 

Different types of crystalline lattices are described by different sets of lattice vectors. 
For instance, one lattice, with the set {a} can be simple cubic, with the particle density 
pi(r + a) = pi(r), while another lattice, with the set {b}, can be face-centered cubic, having 
the density periodicity p2(r + b) = p2(r). 

Gaseous and liquid phases both are characterized by uniform density. Hence the latter 
cannot be an order parameter for distinguishing between these phases. The difference be- 
tween liquid and gas can be noticed on the microscopic level by studying the pair correlation 
function displaying a kind of short-range order in the case of liquid. This type of short-range 
order can be captured by means of the order indices [25||26]. 

Systems with Bose-Einstein condensate and without it can be distinguished by the related 
condensate density po. When po > 0, the system is Bose-condensed, and when po = 0, it is 
not Bose-condensed. 

Superfluid and normal fluid are distinguished by superfluid density ps- The system is 
superfluid if ps > 0, while it is normal if ps = 0. 



Dielectric and metallic phases can be distinguished by the single-electron energy gap AE, 
with AE > representing dielectric phase and AE = corresponding to metallic phase. 

Superconductor and normal conductor are characterized by the presence or absence of a 
gap A in the spectrum of collective excitations. For a superconductor, A > 0, while for a 
conductor, A = 0. 

A system can be described by several order parameters corresponding to different prop- 
erties. For example, it can be ferromagntic and dielectric or ferromagnetic and metallic. 

Thermodynamic phases, having different order parameters (or order indices), are usually 
related to different types of symmetry. Under symmetry one understands invariance under 
some transformations forming a group. A symmetry group is defined as follows. 

Let Qi be a transformation labelled by an index i. And let the family of such transfor- 
mations satisfy the additivity property Qi ■ Qj = gi^j. If the family includes a unit element 
% = 1 and, for each g^ there exists an inverse element, such that g~^gi = 1, then the family 
{c/i} forms a group. 

As an example, consider the inversion symmetry that can be illustrated by inverting the 
magnetization vector M of a ferromagnet. Denoting the inversion transformation by /, we 
have /M = — M. If the order parameter is given by | — M|, it is symmetric with respect to 
inversion: IM = \ — M| = |M| = M. Two elements, the inversion transformation and unity, 
form the inversion group {/, 1}. 

Another example is the rotation symmetry. Denote the spatial rotation transformation 
by R. The order parameter M is invariant under the rotation group {R}, since RM = 
\RM\ = |M| = M. 

One more example is the translation symmetry. Let {T(r)} be a translation transforma- 
tion by a vector r. If a is a vector of a crystalline lattice, then the crystal density is invariant 
under the transformation 

f (a)p(r) = p(r + a) = p(r) . 

The family {T(a)} composes the lattice translation group. 

In the case of a liquid, one has T{r)p = p for any vector r. Then the family {T(r)} forms 
the total translation group. 

2.4 Phase transitions 

A phase transition is a transformation of one thermodynamic phase into another. The change 
happens under a varying thermodynamic parameter, say, temperature T. When T reaches 
the stability boundary of one phase, the system becomes unstable and, under the following 
variation of T, transforms into a stable phase. The point of a phase transition is the point 
of instability. Therefore, at this point, some susceptibilities may diverge. Since thermo- 
dynamic phases can usually be characterized by order parameters, at the phase transition 
point, the order parameter qualitatively changes. Depending on the kind if this change, one 
distinguishes three types of phase transitions. 

Let an order parameter be denoted by r], being in particular cases either magnetization 
M, or particle density p, or condensate density po, or superfluid density ps, or single-particle 
gap AE, or the gap in the spectrum of collective excitations A. And let temperature T be 
varying, with To being the point of a phase transition. The latter are classified into three 
types. 



(i) First-order phase transition. This is a discontinuous transition, with a finite jump of 
an order parameter rj at the transition point: 

\r]{To + 0)-viTo-0)\>0. 

The values of ?7(To ±0) are not necessarily zero. An example of a first-order phase transition 
is the melting-crystallization transition between a crystal and a liquid, where there occurs a 
discontinuity of the density. 

(ii) Second-order phase transition. This is a continuous transition, with a continuous 
variation of an order parameter between its non-zero and zero values: 

\v{n + 0)-v{Tc-0)\ = 0, 

where r){T) becomes zero for one of the phases, say for T > T^, being non-zero for the other 
phase. The point of a second-order phase transition is also called the critical point, since 
some of the susceptibilities diverge there. Respectively, the transition temperature is denoted 
as Tc, to distinguish it from the first-order transition point Tq. Examples of second-order 
transitions are different magnetic transitions or Bose-Einstein condensation. 

(iii) Crossover transition. This is a fast but gradual variation of an order parameter that 
is non-zero from both sides of the crossover point Tc. The crossover point can be defined as 
the inflection point of the order parameter: 

|l = (T = Tc). 

That is, the first derivative of the order parameter displays an extremum at the infiection 
point. 

In the vicinity of a phase-transition point, the system properties may essentially change, 
because of which a region around a transition point, where such changes are strong, is termed 
the critical region. 

If thermodynamic phases are described by different symmetries, the latter also change at 
the phase transition point. When the symmetry of one of the phases is lower than that of 
the other phase, the transition to the lower-symmetry phase is called spontaneous symmetry 
breaking. It is termed spontaneous, since the symmetry breaking occurs in a self-organized 
way, without action of external fields. The cause for the spontaneous symmetry breaking 
is that the phase with a lower symmetry happens to be stable, when the symmetric phase 
becomes unstable. 

More information on phase transitions can be found, e.g., in the books P7H52] . 

2.5 Heterophase mixtures 

In addition to pure thermodynamic phases, there can arise mixtures of different phases, 
where the latter are not spatially completely separated, but are intermixed in such as way 
that each of the phases occupies randomly distributed in space regions of mesoscopic size. 
Such mesoscopic mixtures are also called heterophase mixtures [33]. For example, we can 
imagine paramagnetic regions in a ferromagnet, ferromagnetic clusters in a paramagnet, 
disordered regions in a crystal, crystalline clusters in a liquid, gas bubbles in a liquid, and 
liquid droplets in a gas. 
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For illustration, let us keep in mind the case of a heterophase mixture of two phases, 
characterized by different order parameters, rji and 772- When the system is composed of two 
macroscopic thermodynamic phases, they are said to be in thermal equilibrium with each 
other if their temperatures coincide, Ti = T2. The phases are in mechanical equilibrium, if 
their pressures are the same. Pi = P2. And they are in chemical equilibrium, if their chemical 
potentials are equal, /ii = /i2- A mesoscopic heterophase mixture, strictly speaking, is quasi- 
equilibrium, because of which the above equilibrium conditions are valid only on average [33] . 

Let Vu be the volume occupied by a i^-phase and N^ be the number of particles in 
the i/-phase. Mesoscopic phases can be separated by means of an equimolecular separating 
surface [3ll[35]. In that case, all extensive quantities are additive, e.g.. 



The geometric weights of each phase are quantified by 

u'u = y {iy = l,2). 
This quantity has the meaning of geometric phase probability and satisfies the conditions 

< w,^ < 1 , wi + W2 = I . 

Equations for the phase probabilities are obtained from the minimization of a thermo- 
dynamic potential under the above normalization condition. Defining, for convenience, 

Wi = W , W2 = I — w , 

one has to find the absolute minimum of a thermodynamic potential, say, of the free energy 
F = F{wy. 

F = abs min F{w) . 

w 

The minimization has to include the consideration of w at the boundary, where either w = 
or w = 1. That is, it is necessary to compare three quantities: F{w), with < w < 1, and 
the boundary values -F(O), and F{1). 

The condition of heterophase equilibrium reads as 

dF _ 
dw 

And the condition of heterophase stability is 

d'^F 



dw"^^ 



>0 



These equations define the geometric phase probabilities Wy playing the role of additional 
order parameters. The second derivative d'^F/dw'^ has the meaning of the heterophase sus- 
ceptibility that is required to be non-negative and finite for guaranteeing the system stability, 
similarly to stability conditions (25). 
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2.6 Landau potential 

If the considered phase transition is of second order or of weak first order, when the order 
parameter in the vicinity of the transition point is small, it is possible, assuming that the 
free energy is an analytical function of the order parameter, to expand the free energy in 
powers of the latter, as suggested by Landau [36] . 

Let us consider the critical region of temperatures T ^ Tc, where the order parameter 
7] = 1]{T) is small. The order parameter rj is assumed to be real- valued and such that 
r){T) 7^ for T < Tc, but r]{T) = for T > Tc. The free energy F = F{ri) is expanded in 
powers of the order parameter as 

F = F(0) + aiT] + a2rf + a^rf + a^rj'^ . 

To have a minimum at finite ?7, it is necessary that 04 > 0. 
In equilibrium, one has 

dF 

-p— = ai + 2a2'r] + Sa^r]'^ + ia^r]^ = , 

OT] 

which requires that ai = 0, when rj = 0. The condition of stability is 

d^F I 2^ 

2 [0,2 + 3a3?7 + 60477 j > . 



drf- 

There can happen two cases, with and without inversion symmetry. 

In the case of inversion symmetry, the free energy is such that F{—ri) = F{r]). Then 
03 = and the condition of equilibrium reads as 

T]{a2 + 2aArj^) = . 

From here, we see that the order parameter at low temperatures can take one of the two 
values 



And above Tc, one has 






77 = {T>Tc). 


The stability condition yields 






02 + 60477^ > , 


which implies that 






a2(T) < {T <T, 




a2{T) > {T>T, 



This suggests that one can accept the form 

a2{T) = a{T-Tc). 

The free energy is symmetric with respect to the inversion of rj to —rj. The two order 
parameters have the same magnitude, |r7+| = |?7_|. This situation corresponds to macroscopic 
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degeneracy. When temperature is lowered, the system passes from the symmetric phase, with 
77 = 0, to the phase with a non-zero value of the order parameter, taking spontaneously either 
r/+ or ?7_. This is an example of spontaneous symmetry breaking. 

When there is no inversion symmetry, then 03 7^ 0, and the condition of equilibrium 
defines the order parameter as 



—803 ± a/Qoo — 320204 
V± = ^ • 

To have a real-valued order parameter, it is necessary that 

9al — 320204 > . 
This situation corresponds to a first-order phase transition. 

2.7 Ginzburg-Landau functional 

In the Landau expansion, the order parameter is treated as a scalar. In general, it can be an 
order function rj = //(r, t) of the real-space vector r and time t. The generalization to this 
case has been proposed by Ginzburg and Landau [37] . 

One defines an effective Hamiltonian density H{r) = H{ri{r,t),T, . . .) depending on 
the order function and thermodynamic variables, such as temperature T. The effective 
Hamiltonian 

H = I H(r)dr 



depends on T, reminding rather free energy. This effective Hamiltonian is the Ginzburg- 
Landau functional of the order function. 

The dependence of the order function on the spatial and temporal variables implies that 
the system is not uniform and, in general, not equilibrium. The evolution equation for the 
order function is 

Stationary states can be represented as 

r/(r,t) = ^(r)e-^^*, 

where e is a real-valued quantity. Then the evolution equation reduces to the stationary 
form 

^^^^ 

dip* 

Introducing the grand effective Hamiltonian, or grand potential 



Q = H-e I \^\^ dr 



the stationary equation can be rewritten as 



6ip* 
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which is a condition of equihbrium. 

The function ip can be normahzed to the number of particles A^ by the normahzation 
condition 



l\i:\'dr = N, N = -^. 



This shows that e plays the role of a chemical potential. 

The dynamic stability of the stationary solution is checked by studying small deviations 
from this solution. For this purpose, the order function is written as 

r/(r,t) = [V>(r) + 5^(r,t)]e-*, 

with the deviation 

representing the wave function of collective excitations. Substituting the above representa- 
tion into the evolution equation yields an equation for u. The condition of dynamic stability 
requires that the admissible u be positive (non-negative). This is because ii u = ui + iu2 is 
complex-valued, then 

-ia;t _ -iuiit u)2t iuit _ iuiit -ui2t 

and there appears a divergence for any 002 7^ 0, as time increases, which means that the 
considered solution is not dynamically stable. 

2.8 Superconductor model 



Ginzburg and Landau [37j illustrated their approach by a model of superconductor. In the 
effective Hamiltonian H = J H{r) dr, the Hamiltonian density is written as 

where m is the effective mass equal to two electron masses, m = 2mo, and the effective 
charge equals two electron charges, e = — 2eo; momentum operator p = — iV, and A is 
the vector potential. The quantities U and $0 can be functions of thermodynamic variables 
T,P, . . .. The order parameter rj is such that |?7p gives the superfluid electron density. 
In the presence of an external magnetic field Bq, the effective Hamiltonian becomes 



H 



■^« - ^ 



dr 



The vector potential defines the magnetic field B = V x A and satisfies the Coulomb 
calibration V ■ A = 0. In equilibrium, one has 

6H 

. 



5r]*{r) 
This leads to the equation for the order function 



1 / e \ 2 

;— p A) ri + Uri + %\ri\'^r] = . (27) 

2m \ c / 
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The extremization with respect to the vector potential, 

'" 0. 



can be rewritten as the Euler-Lagrange equations 

6H\ _sr^ d f 6H 

in which A'^a = dA°'/dr'^. From here, one finds the equation 

c 



with the supercurrent 



2 
ifc- / j,_ . . ^v o I |2 A 

ir] Vrj — rjVri ) — — ml A . 

2m mc 



The equation for the order function r/ is a kind of the nonhnear Schrodinger equation. 
It is worth emphasizing that the order-function equation, being a kind of the nonhnear 
Schrodinger equation, has appeared, first, in the paper by Ginzburg and Landau [37J. There- 
fore such equations can be called the Ginzburg- Landau order-function equations. 

2.9 Superfluid vortex 

For systems of neutral Bose atoms, such as superfiuid helium, the effective Hamiltonian 
density takes the form 

where m is atomic mass and the quantities U and $o can be functions of thermodynamic 
variables T, . . .. The order parameter, or order function t] = T]{r,t) gives the density jr^p of 
atoms in the ordered state. 

With the nonlinear Hamiltonian 

H[r]]^-^ + U + <^o\v\\ (28) 

the evolution equation for the order function 

^ ^ = H[v]v (29) 

is the nonlinear Schrodinger equation, being a particular case of the Ginzburg-Landau equa- 
tion (27). For an equilibrium system, one gets the stationary equation 

H[ip]ip = eip . (30) 



14 



The latter equation possesses vortex solutions [38]. To show this, one passes to cylindrical 
coordinates r = {r±,(p,z}, where r± = a/x^ + y'^ is the transverse spatial variable. In this 
representation, 

2 92 Id 1 92 /^2 



dr\ r_L dr± r'j_ dip'^ ' dz"^ 

Looking for the solution in the form 'ip{r) = V'±('"±)c*'^'^, with 'ip±{rj_) real- valued, and im- 
posing the periodicity constraint 

i>{r±, (f + 27r, z) = V'(r_L, (p, z) , 

one finds that the quantum number z/, called the winding number, has to be quantized, 
i/ = 0,±l,±2,.... 

The stationary equation (30) yields 

d? I d v^ 

" ' ''^ V^^ - 2m$oVi = . (31) 



dr\ r_L (ir_L r\ 

Considering a cylindrical system of radius i?, length L, and volume V = irR'^L, the solution 
to this equation is normalized as 

27rL / \ilj^{r^)\\^dr^ = N . (32) 

Jo 

One introduces the healing length 

e = , -^ . (33) 

^Imie - U) 

Using this scale, one defines the dimensionless radius and the dimensionless order function. 



^"e 



""^ f{r)=^^/2m^o^±{r±). (34) 



Then Eq. (31) reduces to 



d"^ 1 d y^ \ „ „-x , ■. 

5;. + ;^ + i-7.)/-/' = 0- (35) 



And the normalization condition (32) becomes 



m$ 



JO 



f\r)rdr = N. (36) 



Looking for a nontrivial solution / 7^ 0, we see that it depends on the winding number v. 

If z/ = 0, there is no vortex, which corresponds to the ground state, when /o = const. 
From Eq. (35), we have /q = 1, hence /o = ±1. This degeneracy is not important, since the 
physical quantity |/oP does not depend on the sign of /q. Without the loss of generality, we 
may choose /o = 1. Normalization (36) gives 

^ = ^T=^' ^ = ^ + ^"^0 (^ = 0). (37) 
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The first vortex solution appears for u = 1, describing the basic vortex. Then Eq. (35) 
can be solved numerically. Close to the vortex axis, where r ^ 1, that is r± <^ ^, the 
solution has the asymptotic expansion 

/(r) ~ VC r(l + C2r^ + C4r^ + Cer'^) , 

with C being a constant and 

_ _ 1 _ 1 + 8C _ _ 1 + 80C 

"'^ ~ ~ 8 ' ''^ ~ ^92~ ' ""' ~ ~ 9216 ■ 

At large distance from the vortex axis, where r ^ 1, that is r± ^ ^, the solution has the 
form 

,. ^ T 1-2 9 _4 161 _6 

•"^ ' - 2 8 16 

Solutions / = fu{r), representing the vortices with higher winding numbers, satisfy the 
equation 

dr"^ r dr r"^ j ^ 

The normalization condition (36) defines the related healing length ^y that is different from 
(37). Notation (32) gives the spectrum 

labelled by the winding number v = 1,2,.... 

If the system radius is very large, then the dependence of the healing length and, re- 
spectively, of the spectrum on the winding number is weak, since, in that case, the main 
contribution to the normalization integral (36) is due to the large region of r ^ 1, where 
/i,(r) — ;■ 1. However for an elongated system, with a small radius, spectrum (38) can essen- 
tially depend on the winding number, if the vortex axis is along the cylinder axis. 

The solution '?/'(r) = '?/'±(r_|_)e*'^'^, for z/ > 0, represents a vortex, since it defines a nonzero 
tangential velocity 

_ J ^ 



e^ , (39) 

p mr 



where 



,2 . ^ 



p^l^l^ j = __(^*v^-^v7A* 



The vortex circulation corresponds to the integral 

V ■ dl = 27r -^ . (40) 

m 

It is worth emphasizing again that the evolution equations for order functions naturally 
appear in the Ginzburg-Landau approach in a similar way for both superconducting and 
superfiuid systems 
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2.10 Landau criterion 

In order to find out whether the system can exhibit superfiuid motion, one employs a simple 
criterion, suggested by Landau [36] . 

Suppose that a flow moves with superfluid velocity v. In the moving frame, where the 
liquid is at rest, the energy per particle is 

E = U + f, 

2m 

with p = |p|. Its energy in the laboratory frame reads as 

with f = |v|. If the energy of an elementary excitation in the moving frame is ^(p), then its 
energy in the laboratory frame is e„(p) = e(p) + p ■ v. 

The fluid motion is dynamically stable when the energy of the elementary excitations is 
positive, 

£^(p)>0, e(p) + p-v>0, 

hence, when — p ■ v < £:(p) for all p. Maximizing the left-hand side of the latter inequality 
with respect to the direction of p, one gets pv < eijp)- Therefore, the superfluid velocity 
should be such that 

V < mm . 

p p 

That is, it has to be lower than the critical velocity 

_ . ^(P) 

Vc = mm . 

p p 

For V < Vc, the superfluid motion is stable, while if v exceeds Vc, the superfluidity is destroyed. 

Strictly speaking, superfluid flow is not absolutely stable, but metastable. The absolutely 
stable state corresponds to the system at rest. 

As an example, when there is no superfluidity, it is possible to mention the ideal Bose 
gas with the spectrum e{p) = p'^/2m, for which minpe(p)/p = 0, hence Vc = 0, and there 
can be no superfluid motion in the ideal gas, even as a metastable state. 



3 Field representation 



Now the basic mathematical notions will be given, providing the fundamentals for the mi- 
croscopic theory of quantum systems. The flrst such a notion is the so-called second quan- 
tization. Below, we shall describe it by employing the fleld representation, following mainly 
Berezin [Ul and Schweber H2]. 
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3.1 Fock space 

Let us enumerate the number of particles by n = 0, 1, 2, . . .. Each particle is characterized 
by a set x of variables from the variable space {x}. For the general mathematical exposition, 
the nature of the particle variables is not important. In particular, this can be the real-space 
vector r, or momentum vector p, or the combination of the real-space vector and a spin 
{r, s}, etc. 

The system microscopic states are the vectors of a Hilbert space. A Hilbert space "H 
is a complete vector space with a scalar product. In this space "H = {/}, for any two 
vectors /, /' G "H there exists a scalar product (/, /'). The latter generates the vector norm 



— a/(/, /). The space is complete, if for a sequence {/'•'^•* G "H | A; = 1, 2, . . .}, the 
limit lim^^oo Z'''^'* = f''°°^ ^ ^- The convergence is understood as the convergence by norm: 
lim,^^||/W|| = ||/(-)||. 

Let us consider the space Tin = {/«}, which is an n-particle Hilbert space composed 
of functions /„ = fn{xi,X2, ■ ■ ■ ,a:„). Let Pij be a permutation operator interchanging the 
indices i and j, with i ^ j, so that its action is defined by the equality 

-* ij Jn \' ■ ■ Xi . . . Xj ■ ■ ■ ) Jn \' ■ ■ -^j ■ ■ ■ -^i ■ ■ ■ ) ■ 

It is assumed that all particles are classified into two types, bosons and fermions. The first 
type is described by symmetric states and the second, by antisymmetric. An n - particle 
function /„ is symmetrized provided that 

In what follows, the upper sign corresponds to bosons, while the lower sign, to fermions. 
The scalar product in Tin is given by the integral 

[jni In) ^ I Jnl-^1) • • • 5 •^n)Jn\-^ii ■ ■ ■ i •^n) dXi . . . uXn ■ 

The integration for each Xi G {xi} is over all {xi}, so that j dx = j\ydx. Each space Tin 

is complete by norm ||/n|| = \/(/n, /«)• In the case of discrete variables, the integration is 
to be understood as summation. For instance, if a; — )■ r, s, where s is a spin index, then 
J dx = ^^ j dr. In that case, the multi-dimensional delta function is 6{x — x') = 6{r — r')6ss'- 

Consider a set of n-particle spaces 'Ho,'Hi, . . ., where the first space "Ho is composed of 
constants /o = const. 

The Fock space is the direct sum 

oo 

J = M) rin ■ 

n=0 



The vectors of the space T = {(/?} have the form of columns 



^ 



h 



In 



In Jn\-^lj ■ ■ ■ 1 •^n) 



The conjugate of ip is the row (/?+ = [/q /i . . . f*^. . ]. For brevity, we shall denote ip = [/„ 
and (/?+ = [/ra]^. The scalar product in the Fock space is 



(V^,V'') 



v^+v?' 



oo 

E 

n=0 



\Jn) Jn) 



for any (p, ip' G J-". Hence, the Fock space is a Hilbert space. The norm is generated by the 
scalar product, \\(p\\ = ^/{}p^^p). 

Recall that functions /„ are assumed to be symmetrized, which can be done by means of 
the symmetrization operator S that makes of functions either symmetric or antisymmetric 
combinations, depending on whether bosons or fermions are considered. Thus, a function of 
two variables is symmetrized as 

Sh{x.,x,)^\[hix.,xj)±hix„x.)]. 

Similarly, functions of many variables are symmetrized. 

Sometimes, one assumes that the symmetrization is done over the products of the type 
f{xi) . . . f{xn). This, however, is not necessary, and the functions fn{xi, . . . ,Xn) can be of 
any type, provided they are symmetrized. 

The Fock space J-', containing only symmetrized functions /„ is called the ideal Fock 
space. 

The vacuum in the Fock space is the column 

1 





10 > 



The restricted Fock space J^n for the fixed number of particles A^ consists of v^at G J^n 
such that (pN = [SnNfN]- The restricted Fock space is the number-conserving space. This 
restriction makes operations in the Fock space more complicated, because of which one 
prefers to deal with the non-restricted Fock space. In what follows, we shall deal with the 
ideal Fock space that, for short, will be called just the Fock space. 
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3.2 Field operators 

The field operators act on the Fock space. One terms ■?/'''" (x) the creation operator and iplx) 
annihilation operator. The action of the creation operator is given by the formula 



ip'^{x)(p 







foS{xi - x) 

\/2 Sf{xi)6{x2 — x) 



n Sfn~l{xi, ..., Xn-l)S{Xn - x) 



In what follows, the latter form will be represented as 

V^^(a;)v9 = ^/n Sfn^l{Xi, ..., Xn^i)5{Xn - x) 

The annihilation operator acts according to the rule 



^(x)^ 



y2/2(xi,x) 
V^ h{xi,X2,x) 



Vn+l fn+liXi, ...,Xn,x) 



This will be denoted as 

llj{x)ip = y/n+1 fn+l{Xi,...,Xn,x) 

The field operators act on the vectors of the Fock space, but ■?/'''" (x) and iIj{x) are not defined 
in the Fock space, since ip'^{x)ip lies outside J-". In that sense, the field operators are the 
operator distributions. 

In the Fock space, one can define the operators 

■?/'t = / /(x)-?/'^(x) dx , tpf = f*{x)il){x) dx , 



where ||/|| < oo. Their action is given by the expressions 

^\ip= Vn Sfn-l{Xi,...,Xn-l)f{Xn) 
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i)fip= Vn + l fn+i{xi,...,Xn,x)f*{x) dx 

which are in J-". 

The action of the annihilation operator on the vacuum is 

'^{x)\ > = . 

And the creation operator acts on the vacuum as 



ip\x)\ > 





6{xi — x) 




A double action of the creation operators yields 



ij\x[)^\x',)\0> 






\/2 S5{xi — x[)5{x2 — X2) 




Acting on the vacuum by n creation operators gives 



l[^\x[)\0> 



i=l 








^' SlY^=iSixi - x'i) 



Therefore, any ip E J^ can be created by the creation operators by the formula 

00 -. „ n 

f = y] -^= / fn{Xl,---,Xn)Y\'^\^i) d^i\ > • 
n=0 ^ 1=1 
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The double action of the annihilation operators gives 



ij{x)ij{x')ip = A/(n + l)(n + 2) fn+2{Xi, ...,Xn,X, x') 
This should be compared with the double action of the creation operators 

■?/'''" (x)-?/'^(x')(^ = 



^n{n - 1) Sfn-2ixi, . . . ,x„_2)5(a;„_i - x)5{xn - x 

Interchanging the variables x and x', we come to the relations 

■ip{x')i/j{x)ip = ±ip{x)ip{x')ip , 

^'^{x')^'^{x)ip = ±tp''(x)ip'^{x')ip . 
Since these equations are valid for any (/? G J-", we have 

[^{x),^{x')]^=[^\x),^\x')]^ = 0. 

Here the notation for the commutator and anticommutator of two operators is used, [A, B]^ 
AB T BA. 

Similarly, we find 

i){x)i)\x')ip = ±{n + l)Sfn{Xl, ..., Xn-1, x)5{Xn " x') 



From here. 



ij\x')ij{x)ip = nSfn{Xi, ..., Xn-l,x)S{Xn - x') 



ijfijlv = [{n + l)U]{f,g) 



^l^fV= [±^/n](/,t/) . 



This can be rewritten as 



Since these equations are valid for any (/? G J-" and any f,g E "Hi, we come to the commutation 
(anticommutation) relations 



[X 



), 2/j\x')]^ = S{X - X') . 



A vector (/? G J-" is called a quantum state or microstate. The Fock space J-" is a space of 
quantum states or a space of microstates. As is shown, the field operators ip{x) and ip^i^x) 
act on the Fock space and can create any vector [/„(xi, . . . ,x„)] in this space. Therefore, 
one tells that the Fock space is generated by ip and denote this Fock space as J-'{ip). 

From the variables x, one can pass to any other required variables, for instance, changing 
the representation from x to k. For this purpose, one defines a complete orthonormal basis 
{ipk{x)} and expands the field operators over this basis. 



^(^) = Yl ^kVk{x) , ip\x) = Y^ al<^l\ 



X] 



The operators a^ and al generate the Fock space J-'{ak). The vectors of the latter space have 
the form [hn{ki, . . . , kn)] G J^{ak). 
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3.3 Operators of observables 

In quantum mechanics, observable quantities are represented by operators of observables. 
An n-particle operator is written as 



n 1 " 

i=l i^j 



i ' j J ~r • • • • 



It is a Hermitian or self-adjoint operator, such that A^ = An. The kernels An{xi, . . . ,a;„) 
are quantum-mechanical operators that are also self-adjoint. Examples of the quantum- 
mechanical operators are the operator of momentum 



Pn = 5Z(-^V.) , 



j=l 



where Vj = d/dvi, the kinetic-energy operator 

2m J 



*»-E(-S). 



and the potential-energy operator 



n 1 " 

in which the pair interaction potential is symmetric, ^{xi,Xj) = $(xj,Xj). The operator of 
the total energy is the quantum- mechanical Hamiltonian H^ = Kn + Wn- 

In the field representation, the general form of an operator acting on the Fock space is 

°^ 1 /" 

A = ^ — ij^Xn) ■ ■ ■ ^''{Xl) An{Xi, ...,Xn) 1p{Xl) . . . ^(x„) dXi . . . dXn , 
^ — n * ^ 



n=0 



with the kernel An{xi, . . . ,Xn) being a quantum-mechanical operator. The operators of 
observables are self-adjoint, A"*" = A, since A'^ = An- 

Examples of the operators of observables in the field representation are the number-of- 
particle operator 

N = ^jj^{x)tlj{x) dx , 



the momentum operator 

the kinetic-energy operator 

K = J^p\x)(-^^^{x)dx, 
and the potential-energy operator 

W = ^Ij\x)U{x)^Ij{x) dx + - / tp'^ {x)tlj\x')^{x, x')tlj{x')tlj{x) dxdx' 
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The total-energy operator is the Hamiltonian H = K + W that reads as 



H= ^^(x) 



2 



^ + ^w 



tlj{x) dx + - I tlj\x)^lj\x')^{x, x')tlj{x')tlj{x) dxdx' . 



The operators of observables act on the Fock space as 

N^ = [nfn] , P^ = [Pnfn] , 
K^ = [Knfn] , W^ = [Wnfn] , H^ = [HJn] ■ 

In the resuhing vector columns, one has quantum-mechanical expressions for n = 0, 1, 2, . . .. 

3.4 Statistical operator 

The system statistical state is characterized by a statistical operator p(t), which is a function 
of time t > 0. This is a self-adjoint operator, p^(t) = p(t). It is semipositive and bounded, 
in the sense of the inequalities 

for all (^ G J-" and v^ 7^ 0. It is normalized to one, Trp(t) = 1, with the trace over the Fock 
space. 

The trace can be taken over any basis {(pk) in ^i which is orthonormal, ^~l^p = Skp, and 
complete, YlkVk'Pk ~ i? where 1 is the identity operator in J-". The trace of an operator A 
is given by TiA = Y,k ft^Vk- 

Observable quantities are defined as the statistical averages 

(i(t)) = Trp(t)i . (41) 

One considers finite systems, with a finite support in real space. Therefore, the operators 
of observables A are termed the operators of local observables. The algebra of local observ- 
ables A = {A} is the family of all operators of local observables describing the treated system. 
And one calls the statistical state the set of all observable quantities <^>={<A>}. 

Being defined for a finite system, with r in volume V, the operator A could be written 
as 

Ay = V" — / ip^Tn) ■ ■ ■ V^^(ri) A„(ri, . . . , r„) V^(ri) . . . V^(r„) dvi . . . dvn , 

where the integration is over the finite system volume. Respectively, the observable quantities 
are < A >y. The average number of particles in the system is A^ =< A^ >y. 
When taking the thermodynamic limit 

A^ 
N ^ 00 , V— j-oo, — — 7- const 

for extensive quantities that are proportional to A^, one should consider the ratio 



lim — < A >v 
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< 00 . 



In what follows, the index V in the integrals and averages is omitted, though always keeping 
in mind that the integration is over the system volume V. 
One tells that two operators 

A= /"i(r)rfr , B= f B{r)dr , 

with the densities A(r) and B{r), satisfy the clustering property if 

< A{r)B{r') > ~ < i(r) >< B{r) > (42) 

for |r — r'l — )■ oo. This property is assumed to hold for the densities of any two operators of 
observables. 

3.5 Evolution operator 

In the Schrodinger picture of quantum mechanics, n particles are described by the wave 
function ^„(t) = "^ni^i, . . . ,Xn,t), where t > is time. Particles are termed identical, 
provided that 

|4^n(t)P=|^„(t)|^ 

where Pij is the operator of particle exchange. The symmetry principle 

A,^„(t) = ±^„(t) 

is a sufficient condition for the particles to be identical. The upper sign corresponds to the 
Bose-Einstein and the lower, to the Fermi-Dirac statistics. 
The scalar product is given by 



(^„, K) = I Ki^u • • • , a^n, t)K{xu 



I JUfi^ '^ ) HiJU\ ... LLJUfi 



with the norm generated by this scalar product, ||\E'„(t)|| = ^(\E'„, \1'„). 

The number of particles n is conserved in time, which imposes the condition ||\l'„(t)|| = 
||^„(0)||. The wave function pertains to a Hilbert space "Hn and satisfies the Schrodinger 
equation 

at 
In the Fock space the wave function can be represented as a collection of the Schrodinger 
wave functions 



*o(t) 



^{t) 



^n(t) 



l^nit)] 
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with ||\l/(t)|| = const. The Fock wave function \l/(t) G J-" satisfies the Schrodinger equation 



where 



H<i!(t) = |ff„*„(«)l 



i I «„(*) 



which is equivalent to the set of the Schrodinger equations in ?/„. 

The time evolution of the wave function in the Fock space can be represented by means 
of an evolution operator U{t) defined so that 

^(t) = f/(t)^(0) , ^+(t) = ^+(0)f/+(t) . 

The evolution operator is additive, U(t + t') = U{t)U{t'). There exists the identity f/(0) = 1 
and the inverse U~^{t). From the norm conservation ||\l/(t)|p = ||\l/(0)|p, we have 

^+(t)^(t) = ^+(0)^(0) = ^+(0)f/+(t)f/(t)^(0) , 

which shows that U^{t)U{t) = 1. Hence the evolution operator is a unitary operator, such 
that U^{t) = U^^{t), with U~^{t)U{t) = 1. By these properties, the evolution operators, 
labeled by the time index, form the evolution group {U(t)}. 

The Schrodinger equation for the wave function \l/(t) is equivalent to the Schrodinger 
equation for U(t), 



tj^U{t) = HU{t), ^lu^{t) = -U+{t)H, 



(43) 



where the fact that the Hamiltonian is self-adjoint, H^ = H, is used. 

In general, the evolution operator cannot be represented in explicit form as a function of 
the Hamiltonian. This is possible, only when the Lyappo-Danilevsky condition is valid: 



H. 



H dt 



. 



(44) 



If this is so, then the formal solution for the evolution operator is 



U{t) =exp\-i I H dt 



This means that the commutator 



H, U{t) 







has to be true. 

In particular, for H not explicitly depending on t. 



U{t) = exp l-iHtj . 
But, in general, the evolution operator cannot be represented in such explicit forms. 
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3.6 Liouville equation 

An n-particle wave function in the Fock space can be written as 







Vnit) 









A pure-state statistical operator is given by the product p„(t) = (/9„(t)(/9^(t). It is the 
normahzed operator, with Trp„(t) = 1, which follows from the wave- function normalizations 
||\l'„(t)|| = 1 and ||(/9„(t)|| = 1. The pure-state operator is idempotent, p^(t) = pn(t). 

A mixed-state statistical operator can be written as 

f){t) = ^u;„v'„(t)v9+(t) . 

n 

Since the statistical operator has to be semipositive and normalized, Trp(t) = 1, it follows 
that the weights Wn are to be such that 

<Wn<l , y^ W« = 1 • 

n 

The time evolution reads as 

p(t) = f/(t)5^ti;„V^„(0)v9+(0)f/+(t). 

n 

Denoting the initial statistical operator 

p(0) = 5^ti;„V^„(0)v9+(0), 



we come to the form 



m = u{t)mu^{t) 



defining the time evolution of the statistical operator. The latter form leads to the Liouville 
equation 

I I- p(t) = \h, Pit)] , (45) 



dt 



in which [i, B]=AB- BA. 
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3.7 Heisenberg picture 

The observable quantities are given by the averages < A{t) >= Trp(t)A, which, using the 
evolution operator, can be written as 

Trp(t)i = Trp(0)f/+(t)if/(t) . 

The Heisenberg representation for the operators of observables is 

i(t) = U+{t)A{0)U{t) , (46) 

with the initial condition A{0) = A. Therefore, we have 

Trp(0)f/+(t)if/(t) = Trp(0)i(t) . 

Then the observable quantities can be written in two equivalent forms: 

< A{t) >= Trp(t)i = Trp(0)i(t) . (47) 

The relation (46) can also be written as the Heisenberg equation of motion 

We should notice the difference of the Hamiltonian form in the Heisenberg equation from 
that in the Liouville equation. In the latter, we have H, while the former contains 

H{t) = U+{t)HU{t) . 

Generally, H{t) does not equal H. This happens only when H commutes with the evolution 
operator U, which holds true when the Hamiltonian H satisfies the Layppo-Danilevsky 
condition (44). 

The Heisenberg equation (48) leads to the time evolution of the observable quantities 



t ^ < A{t) > = < \A{t), Hit) 

at L 



> 



In this way, one can use two representations. In the Schrodinger picture, one deals with 
^(i), p(t), and A(0). And in the Heisenberg picture, one employs ^l'(O), p(0), and A{t). 
The generic form of the operators of observables is 

i = ^ — / i)\Xn) . . .i)\xi)An{xi, . . . ,Xn)i'{Xl) ...^{Xn) dxi . . . dXn ■ 
n=0 ' 

The time evolution of the field operator reads as 

^(x, t) = U+{t)ilj{x, 0)f/(t) , (49) 

with the initial condition ip{x, 0) = iplx). This is equivalent to the equation of motion 

t-ij{x,t)= \ij{x,t), H{t)] , (50) 



dt 
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where H{t) = U+{t)HU{t). 

For the Heisenberg representation of an operator of observable A{t) = lf^{t)AU{t), we 
find 

A{t) = ^ — i^\Xn, t)... ij\xi, t)An{Xi, ..., Xn)lp{Xi,t) . . . 1p{Xn, t) dXi . . . dXn ■ 
n=0 ■ "^ 

If the kernel An{xi, . . . ,x„,t) explicitly depends on time, then the Heisenberg equation of 

motion reads as 

d r 1 3A 

z - Ait) = A{t), Hit) + U^t) ^ — Uit) . (51) 



dt 



dt 



Again, one has to pay attention that, contrary to the Liouville equation, here we have Hit). 
The commutation relations for the field operators from Sec. 2.2, can be represented for 
the time- dependent field operators (49). To this end, we see that 

[^(x,t), i^ix\t)]^ = f/+(t)[^(x), Hx')hUit) , 

and similarly for other relations. Therefore, we obtain the equal-time commutation relations 

[7/>(x,t), ^ix',t)]^ = , [^Hx,t), ^Hx',t)]^ = , [V'(x,t), tlj\x',t)]^ = 6ix-x') . (52) 

According to Sec. 2.3, the Hamiltonian reads as 

V2 



Hit)= Hj^ix,t) 



2m 



+ Uix,t) 



Tpix.t) dx+ 



1 



-F - / il)\x,t)il)\x' ,t)(^ix,x')il)ix' ,t)il)ix,t) dxdx' , 



(53) 



so that the equation of motion (50) becomes 






-^^U(.,t) 



^(x,t)+ ^ix,x')ilj^ix',t)iljix',t)iljix,t) dx . (54) 



It is important to stress that the Heisenberg equation of motion (50) for the field operator 
is equivalent to the variational equation 



. 9 ,, , SHit) 
The proof of their equivalence is given in Ref . |13] . 



(55) 



3.8 Local dynamics 

Using the evolution equations in the Heisenberg representation, it is straightforward to derive 
the equations for different local density operators [H]. Let us consider the local density of 
particles 

n(x, t) = '?/'^(x, t)'?/'(x, t) . (56) 
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The density commutes with itself, [n{x,t),n{x',t)] = 0, for any statistics. The integral of 
the density gives the operator of the total number of particles 



N 



h{x, t) dx . 



And let us define the local current density of particles j = i{x,t) by the equation 

[^t(v^) - iW)^] , (57) 



J 



2m 



where, for short, we write ip = ip{x,t) and a: = {r, . . .}. 
First, it is easy to derive the continuity equation 



dfi 
'dt 



+ V-j = 0. 



(58) 



The equation for the current density follows: 

f + i^ {(^V)(V^) - ^tv(V2^) + (V^t)(v2^) _ [V(VV)] ^} = (59) 

1 
where 



771 



V(f/ + <l>) 



^, 



1) = ^(x,t) = / ^{x,x') n{x',t) dx' . 

Analogously, we can get the equation for the time variation of the Hamiltonian density 
(density of energy) 



H{x, t) = —- V^Ux, t) ■ V^(x, t) + ijix, t) 
2m 



U{x) + -^x,t) 



ij{x,t) , 



which can also be written in the form 
H{x,t) =tlj\x,t) 



v^ 1 ^ ■ 

_ + Uix) + -^x,t) 



ij{x,t) 



(60) 



(61) 



The total energy Hamiltonian is H = J H{x, t) dx. 

It is worth emphasizing that the variation in space bears some analogy with the variation 
in time. To show this, let us treat x as r. Consider the operator of momentum 

P = m / j(r) dr = tlj\r){—i\/)^p{r) dr . 

We see that the commutator 

[ipir), P] = -iV^(r) 

reminds us the Heisenberg evolution equation (50). Introducing the translation operator 

f (r) = exp(ir ■ P) , (62) 
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we find 



^(r) = f +(r)V'(0)f (r) 



This is similar to the dependence of the field operator on time, as expressed through the 
evolution operator in Eq. (49). That is, the translation operator (62) describes the spatial 
variation in a manner similar to the description of the time evolution by means of the 
evolution operator U(t). 

The translation operator (62) enjoys the properties 

f (0) = i , t+(r)f (r) = i , f-^(r) =f(-r) = exp(-ir-P) , t(r)t(r') = t(r + r') . 

Hence the family {T(r)} forms a unitary translation group. 

3.9 Coherent states 

Coherent states are the eigenstates of the annihilation operator. A coherent state \ri >E J^ 
is a solution of the eigenproblem 



ip{x)\ 1] > = ri{x)\ 1] > 



(63) 



where iIj{x) = iIj{x,0) and ri{x) is called coherent field. Being a vector of the Fock space, the 
coherent state has the form 



r]i{xi) 
^2(a;i,a;2) 



|?7 >= 



^n\Xij • • • 5 Xn) 



For short, we shall denote this coherent state as 

1^ >= ['nn{Xi,...,Xn)] ■ 

A conjugate to \ri > is \ri >+=< ri\. The scalar product of two coherent states is 



<^l^>= ^{Vn,^n) , 



n=0 



where the scalar product of the components is 



yVni <,n) 7 VnK'^ii • • • i •^n)<,n\-^li ■ ■ ■ ■> •^n) CLXi . . . CLX^, 
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By the definition of tlie anniliilation operator, we liave 



ij{x)\r]>= Vn+lr]n+i{xi,...,Xn,x) . 
Taking into account tlie definition of tlie colierent state (63), we get tlie recursion relation 



Vn+l r]n+iixi, ...,Xn,x) = r]{x)r]n{xi, • • • , a;„) . 
Tlie iterative solution of the latter recursion relation results in the form 

n 



r/„ (xi,..., x„) = —^ TT r7(a;i) . (64) 

Vn! fj^ 



The scalar product of the components (64) gives 



n! 



The scalar product of two coherent states reads as 



n=0 

where 

{V,0 = / V*ix)^ix) dx . 

The coherent state can be normalized as 

II \ri > IP =< ri\ri >= 1 , 

which defines 

\t]o\ =exp|- - {r],r]) 

In this way, the scalar product of two normalized coherent states becomes 

r 1 1 

< ^IC > = exp <^ - - {r], T]) + [t], i)- - (^, 

Since 0<|<?7|.^>|<1, it follows that the coherent states are not orthogonal with each 
other. The set {\ri >} forms an overcomplete basis for which the resolution of unity holds 
true: 

/ \ri >< r]\ Vr] = 1 . 

Integration over the coherent fields can be introduced as follows. Let us take an orthonor- 
mal complete basis {Xk{x)} and let us expand the coherent field as 



Vix) = ^^CkXkix) , Cfe = |cfc|e*"*= , 
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with ak real- valued and < |ca;| < 1. Then the integration over the complex- valued coherent 
field is defined with the help of the differential measure 

k 

A coherent state jr^ > is an eigenstate of the annihilation operator, but it is not an 
eigenstate of the creation operator ■il)'^{x). This is seen from the action 

^\x)\ 7]> = y/nSr]n-l{Xl,...,Xn-l) S{Xn - x) 

However, we have the equation 

< r]\ ^\x) = < r]\ r]*(x) . 
In view of Eq. (64), the coherent state can be represented as 

oo / „ \ n 



n=0 

In turn, the latter takes the form 



V> = ^ ^ ( vix)i'\x) dx] I > . 



1^ > = Vo^W \ / 'rii^)^H^) dx>\0> . (65) 



Since the coherent state pertains to the ideal Fock space, its components possess the 
symmetry property 

J^^ij^n\X\^ • • • ) -^raj '^'nn\X\^ • • • ; Xn) 

with respect to the permutation of the variables Xj and Xj, for i ^ j. Therefore the coherent 
fields satisfy the commutation relations 

r]{xi)r]{xj) - 7]{xj)r]{xi) = (Bose) 

in the case of bosons and the anticommutation relations 

rj{xi)rj{xj) + rj{xj)rj{xi) = (Fermi) 

in the case of fermions. For bosons, the coherent field is a complex-valued function. While 
the fermionic coherent fields are Grassmann variables. 

Time-dependent coherent states can be defined in the same way as it is done for other 
states in the Fock space: 

\ij(t) > = U{t)\r]>=[r]nixi,...,Xn,t)] . (66) 

The invariance of the normalization is assumed, so that 

< r]{t)\r]{t) > = <r]\r] >= 1 . 
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By treating the coherent state as an eigenstate of the annihilation operator satisfying the 
eigenproblem 

^PixMt)>=r]{x,tMt)>, (67) 

we get the definition of a time- dependent coherent field ri{x,t). The latter equation is 
equivalent to 

ij{x)U{t)\ri > = ri{x,t)U{t)\ri > , 

from where 

U^{t)^lj{x)U{t)\r] > = r]{x,t)\r] > . 

With ilj{x,t) = U^{t)ilj{x)U{t), we obtain another form of Eq. (67) as 



ip{x,t)\ri > = ri{x,t)\ri > . 



(6^ 



Using the time- dependent coherent states, it is possible to pass to the coherent-state 
representation for the operators of observables. For instance, the number-of-particle operator 
yields 

< r]\N{t)\r] > = \v{x,t)\'^ dx . 



And for the Hamiltonian, we find 

< ri\H{t)\ri > = I rf{x,t) 
1 



— + Uix) 



ri{x,t) dx+ 



+ ^ / \'n{x,t)\'^^{x,x')\i]{x',t)\'^ dxdx' . 
The temporal derivative gives 

, . (9 , , , . ■ 9 , , 

<7]\l — ^(x, t)\r]> = i — 7]{X, t) . 

The evolution equation (55) transforms into the equation 

.8 S <r]\H{t)\7]> 

Explicitly, we obtain the coherent-field equation 






V^ f 

h U(x) + / ^(x,x' 

2m J 



)\ri{x,t)\^ dx' 



r]{x,t) . 



This is the nonlinear Schrodinger equation. 
Presenting the solution in the form 

r/(x,t) = ?7(a;)e"*'^* 

we get the stationary coherent-field equation 



+ U(x) + / $(a;, x') \r](x')\'^ dx' 

2m J 



ri{x) = Eri{x) . 



(69) 



(70) 



(71) 
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It is important to keep in mind that, as is seen, the above equation has a meaning, provided 
that the interaction potential is integrable, such that 



$('x,x') dx' 



< oo . 



The coherent-field is defined by this equation up to an arbitrary global phase, that is, the 
general solution is ria{x,t) = e^°'ri{x,t), with a being an arbitrary real-valued quantity. 
More information on coherent states can be found in the books 

3.10 Density matrices 

Statistical averages of the products of field operators define the reduced density matrices [25] . 
An n-order density matrix is defined as 

p„(xi, ...,Xn,x[,...,x'^,t) = Tr i>{xi) . . . tlj{xn)pit)ij\x'^) . . . i>\x[) = (72) 

= Tr i'{xi,t) . . . ip{xn, t)pip^{x'^, t)... i'^{x[,t) = 

= < tfj^x'^, t)... 'llj\x[,t)'llj{xi,t) . . . ^piXn, t) > . 

They are symmetric with respect to the simultaneous permutation of Xj with Xj and x'^ with 
Also, they enjoy the following properties: 

Pnl-^l' • • • ? •^ni 2^1; • • • 5 ^ni ''/ Pny^ni ■ ■ ■ i-^li -^ni ■ ■ ■ i ^li <^) i 
Pn\Xl) • • • 1 Xni -^l; • • • ; Xrf^i i) PfiyXfiy ■ ■ ■ i Xi, 3;^, • • • , X-^^, Zj , 
PnV^^i • • • 1 •^ni ^li ■ • ■ 1 ^ni ''/ Pny^li ■ ■ ■ i ^ni -^li ■ ■ ■ i -^n-, <^) ■ 

The matrices are called reduced because their partial integration reduces them to lower- 
order matrices. Suppose we consider the restricted number- conserving space J-jv, for which 
the number of particles is fixed, so that for ip^ G J-jv, one has Nip^ = Nip^. Then, it follows 
that < NA >= N < A >. For the second-order density matrix, we have 

P2(a^, x", x', x", t) = < h{x", t)ip\x', t)ip{x, t) > — 5{x' — x")pi{x, x", t) , 

where n{x^t) = il)'^{x^t)ip{x,t) is the number-of-particle density. Integrating over x" yields 

P2{x,x",x',x",t) dx" = < Nip\x',t)ip{x,t) > -pi{x,x',t) . 

On the restricted Fock space J-jv, 

< Ni:^{x',t)ij{x,t) > = Npi{x,x',t) . 

Therefore, 

P2{x,x", x',x",t) dx" = {N — l)pi{x,x' ,t) . 
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The name of the matrices contains the term density, since the diagonal element of the 
first-order matrix 

Pi{x,x',t) = < il)\x' ,t)il){x,t) > 

defines the particle density 

Pi{x, X, t) =< ip'^(x, t)^{x, t) >= p{x, t) . 
The name matrices means that these objects can be treated as matrices. Thus, the matrix 

Pn [Pn\^lj ■ ■ ■ 1 ^ni X^^, . . . , X^, tj\ 

is a matrix with respect to the x-variables. This matrix is evidently self-adjoint, p^ = pn, 
and semipositive, since 

[jn-i Pnjn) 

Here, the scalar product 



>0. 



\Jny PnJn) — / Jn v-^l) • • • 5 •^n)P\-^ly • • • ; -^n; 3)]^, . . . , X^, ZjJnyX-^^y ■ ■ ■ y •^n) '^•^l • • • wX^ UX^ . . . CtX^ 

and the vector norm 

I \In\ I l/n? In) / In [,-^li ■ ■ ■ i •^n) Jx[,-^li ■ ■ ■ i •^n) wXi . . . CLXn 



are used. 

The norm of an n-order density matrix is given by the expression 

llPnII = sup \{fn,pnfn)\ ■ 

ll/n|l=l 

The eigenproblem 

Pnjnj -^nj Jnj ; 

implies the equation 

Pn\Xij . . . , Xnj X-]^, . . . , X^j Zj Jfij\X-^^ • • ■ i X^) uX-^ . . . dX^ ■^nj Jnj\Xi, ■ ■ ■ ■, Xn) . 

Being self-adjoint and semipositive, the density matrix enjoys real non-negative eigenvalues. 
The largest of the eigenvalues defining the norm 

||p„|| = sup A„j = An 
J 

is of special importance characterizing the level of ordering in the system. 
The normalization for the density matrices can be summarized as follows: 



Trpi = / pi(x,x,t) dx = N , 
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Tipr, 



Pn[,-^1^ ■ ■ ■ ! -^ri! -^1; • • • ; -^m '^j wXi . . . CLXn 



N\ 



{N -n)\ ■ 
The level of ordering in the system can be described |25l|26] by the order indices 

log||Pn|| 



(^{Pn) = 



log Trp„ 



(73) 



Keeping in mind that A^ ^ n > 1, we may use the Stirling formula z\ ~ y2Ti e ^ z^"^^'^, in 
which \z\ — )■ oo. Then, we have 



N\ 



{N-n)\ \e 
As a result, the order index (73) reads as 

(^{pn) 



[NY , N\ 

I — 1 , log 



{N -n)\ 



n log A^ . 



log An 

nlogiV 



(74) 



for A^ > ra > 1. 

The order indices allow us to distinguish between long-range order, mid-range order, and 
short-range order [25112^ . This is contrary to the order parameters that distinguish only the 
long-range order from the absence of such an order. 

By means of the reduced density matrices, it is possible to express the statistical averages 
of operators of observables. In the general case, we have 

< A{t)> = y^ — / lim An{xi,...,Xn)pn{xi,---,Xn,x\,...,x'^,t) dxi...dxn. 

In particular cases, for the average of the number-of-particle operator, we get 



< N > = pi{x,x,t) dx = N 



For the kinetic-energy operator. 



V2\ 



< K > = I lim I — - — j pi{x,x',t) dx . 



And for the potential-energy operator 
<W > = 



U{x)pi{x,X,t) dx + - I ^{xi,X2)p2iXl,X2,Xi,X2,t) dxidx2 



The total internal energy is E = <H> = <K> + <W>. 

The equations of motion for the reduced density matrices follow from their definition and 
Eq. (54) for the field operator. The first-order density matrix satisfies the equation 



i -K7 pl{Xi,X2,t) 



Y1^U(..)^^-U(.,_) 



Pl{xi,X2,t) + 



(75) 
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+ / [$(xi,X3) - $(a;2,X3)] p2{xi,X3,X2,X3,t) dX3 , 

involving the second-order matrix. The equation for the second-order matrix involves the 
third-order matrix, and so on. The sequence of these equations is called the Bogolubov - 
Born - Green - Kirkwood - Yvon hierarchy [1^1171 - 39] . 

Particle correlations are characterized by the pair correlation function 

< ^^{x,t)^^{x',t)^{x',t)ij{x,t) > 
^^''''''^" P{x,t)p{x',t) 

which is expressed through the second-order density matrix: 

P2{x,x',X,x',t) 



g{x,x',t) 



p{x,t)p{x',t) 



The pair correlation function is real, g*{x,x',t) = g{x,x',t) and symmetric, g{x',x,t) = 
g{x,x',t). It satisfies the asymptotic condition 

g{x, x', t) — 7- 1 (|x — x'l — 7- oo), 

following from the clustering property (42). 

4 Statistical states 

In this section, different types of statistical states are treated and their properties are dis- 
cussed. Several important mathematical facts are presented, which are widely used for 
describing quantum systems. 

4.1 Gibbs entropy 

Gibbs entropy is the central notion for characterizing statistical states. It is also called 
Shannon entropy or information entropy, since it quantifies the amount of the unavailable 
information on the considered system. For a statistical operator p, the Gibbs entropy is 

5(t) = -Trp(t)lnp(t) . (76) 

Since p{t) is semipositive, bounded, and normalized to unity, S(t) is semipositive, S{t) > 0. 
Here and in what follows, the functions of operators are defined by the rule 

n=0 

where a is a convenient nonoperator quantity. Below, several important facts, related to the 
Gibbs entropy are proved. 

Liouville theorem. Gibbs entropy does not change with time: 

S{t) = SiO) . (77) 
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Proof. The time evolution of the statistical operator is given by the relation 

p(t) = f/(t)pf/+(t), 

where p = p(0). By the definition of the function of an operator, we have 

Inp(t) = U{t) \npU+{t) . 

Therefore 

S{t) = -Trp(0)lnp(0) = S{0) , 

which proves equality (77). D 

Gibbs inequality. Let A and B be semi-positive operators, possessing traces, then 



Trilni-Triln5 > Tr [i-^V (7J 



Proof. For a semi-positive operator A, one has 

lni> 1-A-K 

Similarly, 

In i - In 5 = \n{AB~^) > 1 - A'^B . 

From here, 

iln(ifi-i) >A-B. 

Then inequality (78) follows. D 

Let us introduce the entropy operator 

S{t) = -\np{t) . (79) 

And let us define the average entropy 

'S{t) = Trp(0)^(t) . (80) 

Notice that the average entropy does not equal Trp(t)S'(0). 

Nondecrease of average entropy. The average entropy (80) does not decrease with 
time: 

S{t) > S{0) . (81) 

Proof. By the definition of the average entropy, one has 

'S{t) - :S(0) = Trp [Inp - Inp(t)] . 

The proof follows the same arguments as in the proof of the Gibbs inequality. We have 

Inp — Inp(t) = In pp^^(t) . 

Since 

lnpp-\t)>l-p-'m. 

39 



we get 

Sit) - S{0) > Trp [1 - p-^p{t 

Then 

Trp [1 - r^pit)] = Tip - Trp(t) . 

As far as Trp = Trp(t) = 1, we come to the inequahty 

S{t) - S{0) > , 
which yields Eq. (81). D 

4.2 Information functionals 

The Gibbs entropy serves as an efficient tool for defining the form of the statistical operators. 
This is done employing the principle of minimal information, which implies a conditional 
maximum of entropy under given additional constraints. The idea of this principle goes back 
to Gibbs [34l|35l[50] and has also been considered by Shannon |5T] and Janes |52]. General 
formulation and different applications have been studied in Refs. [331H6] . 

Suppose that some operators of local observables Ai E A at t = are assumed to give 
the known observable quantities 

A = <A> = TipAi , (82) 

enumerated with i = 0,1, 2,.... One of the quantities is always Aq = 1, which guarantees 
the normalization condition 

< i >= Trp = 1 . 

The operators entering constraint (82) are termed the constraint operators. The information 
functional, or simply information, is 

I[p] = -S + Y,Xi (rrpAi - A,>j , (83) 

i 

where A, are Lagrange multipliers, controlling the validity of conditions (82) through the 
variation 6I[p]/SXi = 0. The statistical operator is supposed to be such that to provide the 
minimum of information (83). Minimization with respect to p implies 

5p dj? 

Using the derivatives 

— Trpij = Ai , — Inp = p"^ 
op op 

leads to the Gibbs statistical operator 

p=i expl ->;A,A, 1 , (84) 



i>0 / 
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where 

Z = Trexp 



i>0 / 



is called the partition function. The second derivative is positive, since 

In the same way, it is possible to find the statistical operator, when the local operator 
densities at t = yield the known densities of observables 

A{x) = <Mx)> = Trpi^(a:) • (85) 

Then one defines the local information functional 

Iioc[p] = ~S + J2 [ H^) I^^pMx) - A^x)] dx . (86) 

The minimization of information (86) results in the statistical operator 



— exp <! - ^ / \i{x)Ai{x) dx \ 



p = — exp <( - ^ ^ / \i{x)Ai{x) dx ) , (87) 



with the partition function 



Z = Trexp < — N^ / \i{x)Ai{x) dx > 



The statistical operator (87) is called the local- equilibrium or quasi- equilibrium statistical 
operator. 

4.3 Equilibrium states 

A statistical state is the set {< A{t) >} of the averages of all operators of local observables. 
By this definition, a statistical state possesses several general properties that are studied in 
the algebraic approach [53]. 

An equilibrium state does not depend on time, so that < A{t) >=< A(0) > for all 
averages from the state {< A{t) >}. This implies that pit) = p(0) = p, because of which 
[H, p] = 0. The latter means that the statistical operator p is an integral of motion. Hence 
it has to be a function of H and other integrals of motion. An integral of motion C is an 
operator for which 

[(7, i^] = , ij^C = [C,H]=Q. 

The Hamiltonian H is, clearly, also an integral of motion. 

Replacing in the information functionals Ai by Cj, and minimizing the information, we 
get the Gibbs statistical operator 



1 

— exp 



Y.\c^ 
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with the partition function 

Z = Trexp I — 




The statistical state {< A >}, composed of the averages < A > = Tip A, with the Gibbs 
statistical operator, is termed the Gibbs state. 

One says that p defines a statistical ensemble, which is a collection of possible dynamical 
states of a physical system. More rigorously, a statistical ensemble is a pair {J-', p} of the 
Fock space of microscopic states J-' and a statistical operator p. Depending on the choice 
of the constraint operators, there are different statistical ensembles represented by different 
statistical operators. 

If, in addition to the normalization condition, one considers the sole constraint operator 
that is the Hamiltonian energy operator H, then the principle of minimal information yields 
the Gibbs canonical operator 

p = 1 e-^^ , Z = Tie-^^ , (88) 

where /3 is inverse temperature, (3T = 1. The statistical operator p describes the Gibbs 
canonical ensemble. In this ensemble, the role of a thermodynamic potential is played by 
the free energy F = F{T, V, N, M), which can be defined by the equality 

Tre-^^ = e-^^ . 

Therefore 

F = -T In Z = -T In Tre^^^-^ . (89) 

Other thermodynamic quantities can be found through the free energy, for instance 

^ = -(f),, = ''(<^>-^)' ^-<^>=^ + rs. 

The Gibbs entropy (information entropy) coincides with the thermodynamic entropy, 

S = -TTp\np = I3{E-F) . 

When, in addition to the normalization condition, as the constraint operators one takes 
the Hamiltonian energy H and the number-of-particle operator N, then one gets the grand 
canonical operator 

p = ^^~^(H^,N) ^ Z = TTe-^(^-^^) , (90) 

Zj 

where p is chemical potential. This statistical operator describes the grand canonical ensem- 
ble. The corresponding thermodynamic potential is the grand potential Q = i7(T, V, /i, M) 
defined by the equality 

This means that 

fi = -TlnTre-^(^-'^^). (91) 



42 



The grand potential Q = F — fiN defines other thermodynamic quantities, for example, the 
entropy 

The Gibbs and thermodynamic entropies coincide, so that the above expression of the entropy 
equals 5* = — Trplnp. 

Generally, different statistical ensembles can describe different physical situations. If the 
considered systems are physically different, then one of the ensembles may be appropriate 
for characterizing one of the system, while another ensemble may be not appropriate for the 
same system, but appropriate for the other system. In that sense, the statistical ensembles 
do not need to be completely equivalent. Though in the majority of the cases, the ensembles 
are equivalent, so that any of them can be used. 

It is possible to distinguish the ensemble equivalence in two senses. One meaning is the 
thermodynamic equivalence of ensembles, when the thermodynamic characteristics, calcu- 
lated for different ensembles, coincide, after being expressed through the same variables, so 
that 

Another meaning is the statistical equivalence of ensembles, when the statistical states 
coincide, so that the observable quantities calculated in one ensemble, say in the canonical 
ensemble giving < A >n, coincide with the same observable calculated in another ensemble, 
say in the grand canonical ensemble giving < A >^, after taking into account the relation 
between the chemical potential and the number of particles. 

4.4 Representative ensembles 

In order to correctly characterize a system, one has to chose a representative statistical 
ensemble, which is an ensemble that uniquely describes the considered system [221E11IS11ES] • 
Actually, the problem of ensemble equivalence is meaningful only for representative ensembles 

An equilibrium statistical ensemble is a pair { J^, p} of the Fock space and a statistical 
operator. To uniquely describe the system, the ensemble must be representative, which 
requires that the statistical operator has to be defined taking into account all constraints that 
are necessary for the unique system description. For this, one has to consider the necessary 
constraint operators, C^ = Ci that are assumed to be self-adjoint, but not necessarily 
integrals of motion. These constraint operators define the corresponding statistical conditions 

a =<C^>= TrpA , (92) 

with the trace over J-". 

The first standard constraint is the normalization condition for the statistical operator, 

1 = < i^ > = Trp , 

where Ijf is a unity operator in J-". Another usual condition is the definition of the internal 
energy 

E =<H >= Tr pH . 
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Keeping in mind the normalization condition, the definition of the internal energy, and 
other constraints (92), required for the unique description of the system, we construct the 
information functional 

I[p\ = Trplnp + Ao(Tr/5 - 1) + /3(Trpi^ - E) + /3 J] z/,(TrpA - Q) , 

i 

in which Aq, /3, and /3z/j are the Lagrange multipliers. 

The minimization of the information functional gives the statistical operator 

p = i e-^^ , (93) 

in which the grand Hamiltonian 

H = H + Y,ViCi (94) 

i 

is introduced. The partition function is Z = exp(l + Aq) = Trexp —(5H. Such a statistical 
operator defines the representative statistical ensemble {J-",/)}. 

4.5 Action functional 

A nonequilibrium statistical ensemble is a pair {J-", p(t), where the statistical operator is a 
function of time characterized by an evolution operator and an initial condition p = p(0). 

In quantum field theory, the system evolution is defined by an action functional [56] . The 
temporal energy operator is written in the form 

r Q 

E= 'ip^ix,t) i — tp{x,t) dx . (95) 

The system dynamical properties are described by the Lagrangian L = E — H and the action 
functional 

A[ij] = f (l- UiC,) dt= f (E-H^dt, (96) 

taking into account the constraints that uniquely define the system. Here H is the grand 
Hamiltonian (94). 

The equations of motion are given by the extremum of the action functional, 

MM _^_ IfL^o. (97) 



S^^{x,t) ' dip{x,t) 

From the form of the action functional (96), it follows that the extremization of the action 
functional is identical to the variational equation 

And the latter, as has been proved ^43j, is equivalent to the Heisenberg equation of motion 

i-ijix,t) = [^ljix,t),H]. (99) 
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A representative ensemble takes into account all statistical conditions, necessary for the 
unique description of the considered statistical system. By defining a representative ensem- 
ble, one comes to the notion of the grand Hamiltonian. And, as is shown above, the grand 
Hamiltonian governs both the statistical properties of equilibrium systems as well as the 
temporal evolution of nonequilibrium systems. This means that the Liouville equation for 
the statistical operator must also contain the grand Hamiltonian. 

4.6 Fluctuations and stability 

In section 2.2, it is explained that in equilibrium systems, susceptibilities are to be posi- 
tive and finite. Here we show how the susceptibilities are connected with fluctuations of 
observable quantities. 

Observable quantities are represented by self-adjoint operators. Fluctuations of an ob- 
servable, corresponding to an operator A, are characterized by the variance or dispersion 

var(i) = < i+i > - I < i > |2 . (100) 

Generally, variance (100) can be defined for any A. For a self-adjoint A+ = A, it reduces to 

var(i) = <A'^> - < i >2 . (101) 

By this definition, var(A) > 0, since 

var(i) = < {A- < A>y > . 

That is why susceptibilities are to be non-negative. 
Let us consider specific heat 

^l_fdE\ ^TfdS\ 

^ N{dT)y^ NydTjy^' 

which can also be expressed through the second derivatives 



a 



T fd^F\ T fd^n 



V — ~ TT 



NKdT^Jy^ N\dT\y^ 



Using this definition, for the specific heat we have 



^ var(if) 

Cv = -^ , (102) 



which shows that it characterizes energy fluctuations. 
Consider now the isothermal compressibility 



1 fdV\ 1 fdN\ 



kt = — — 



- 



V \ dP J rp^ Np ydfJ, J j,y 

which can also be written as 

1 fd^n\ 1 fd^P' 

kt = — 



Np \ dp^ ) j.y P^ \ dp^ J rpy 
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Therefore for the compressibihty, we get 



kt = ^}-J , 103 

NpT ' ^ ^ 



which means that it describes the number-of-particle fluctuations. 
The number-of-particle variance can be written as 



var(A^) = [< h{x)h{x) > —p{x)p{x')] dxdx' , 

being expressed through the density- density correlation function 

< h{x)h{x') > = p2{x, x', X, x') + p{x)6{x — x') . 
Therefore, it can be written as 

var(iV) = N + p{x)p{x')[g{x,x') - 1] dxdx' . (104) 
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Considering the magnetic susceptibility, we take into account that an external magnetic 
field B enters the Hamiltonian through the expression — X]j=i ^i ' B, where Mj is the 
magnetic-moment operator of an i-particle. The magnetic susceptibility tensor is defined as 

1 (9M" 

Xal3 



N dB^ ' 
where the components of the average magnetization are 

N 

- = ^ A>° ^ ivr = " 
dB 



M" = - - — = < if" > , M = V Mi 



Then the susceptibility tensor is 

1 d'^G _< M^'M^ > - < M" >< M'^ > 

The diagonal elements of the tensor give the magnetic susceptibilities 

var(M") 
Xaa = ^y ' , (105) 

which shows that they describe magnetic-moment fluctuations. 

Being proportional to the variances of the related observables, all these susceptibilities 
are non- negative: 

Cy > , KT>0, Xaa>0, 

which is a part of the stability conditions. For the stability of an equilibrium system it is 
also required that the susceptibilities be finite: 

Cv < OO , Ky < OO , Xaa < OO . 
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This means that the stabihty conditions, in terms of the variance of observables, are 



varfA) 
< — -V^ < oo. 

- N 



(106) 



This has to be vahd for any A^, including A^ — )■ oo. 

When var(A) oc N, this is called the thermodynamically normal dispersion corresponding 
to thermodynamically normal fluctuations. And if var(A) oc A^''', with 7 > 1, this is called 
the thermodynamically anomalous dispersion describing thermodynamically anomalous fluc- 
tuations. Thermodynamically anomalous fluctuations imply the system instability. 

Summarizing, the susceptibilities 



a 



V 






Kt 



var(iV) 
NpT 



Xc 



var(M") 

NT 



describe the fluctuations of the corresponding observables, they are expressed through the 
variance of the related operators and have to be non-negative and finite for any equilibrium 
system. 

4.7 Gibbs-Bogolubov inequality 

A very important result, widely used in statistical physics is the inequality, first, proved by 
Gibbs for classical systems and then extended by Bogolubov to quantum systems, because 
of which it is named the Gibbs-Bogolubov inequality. 

Let us consider two self-adjoint operators A and 13. And let us define the functional 



F{A) = -T\nTre'^'^ 
We introduce the averages 

<A>A = TtpaA 
in which 

Pa 



F{B) = -T\ni:ie-^^ . 



<A>B = TrpB^ , 



-pA 



-PB 



PB 



Tre-/3S 



Tre-^^ 
Then the following theorem is valid. 

Gibbs-Bogolubov inequality. For the quantities F{A) and F{B), one has 

F{A) <F{B)+ <A-B>B 

and 

<A-B>A + F{B) <F{A) . 

Proof. According to the above definitions, we have 

TTpB\npB = /3\F{B)- <B>B 



(107) 
(108) 



TrpB In Pa = /? F{A)- <A> 



'B 
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By the Gibbs inequality (78), 

TrpB In Pb > Trp^ In Pa ■ 

From here, we get inequahty (107). Interchanging A and B, we come to inequahty (108). D 

The most often used apphcation of the Gibbs-Bogolubov inequahty is for finding the best 
variational approximation for a complicated problem. Suppose we need to calculate the free 
energy 

F{H) = -TlnTre-'^^, 

which is impossible to calculate exactly. Then we can take an approximate free energy 

F(^„pp) = -TlnTre-'^^-% 
defined through an approximate Hamiltonian Happ. We use the notation 

exp{-/3H) 



<A> = TTpA, p 

^ ^ ^app -Lr Papp^ ) 

Then, by the Gibbs-Bogolubov inequality. 



Trexp(-/3/J) 

_ exp{-/3Happ) 
Trexp(— /3i7, 



appj 



F{H) < F{Happ) + <H- Happ >app • (109) 

The best approximation is obtained by minimizing the right-hand side of inequality (109), 
that is, by choosing the approximate Hamiltonian that provides such a minimization. 

A particular form of the Gibbs-Bogolubov inequality is derived for zero temperature, 
when T — )■ 0, hence /3 — ?> oo. By the definition of the trace, 

k 

where (fk can be taken as the eigenfunctions of the Hamiltonian, 

Hipk = Ekfk ■ 

The ground-state energy is 

Eo = min Ek = (v?o, Hipo) . 

k 

Assuming a non-degenerate spectrum, the above trace can be rewritten as 



Eri^-''"^' = 


= e"'^^" 


l + ^e-^iE,-Eo) 


k 




I kytO J 



Since i?fc — -^o > for /c 7^ 0, we have 

Tre"^^ ~ e-'^^" {/3 -^ 00) 
48 



Therefore 

F{H)^E, (T^O). 

Respectively, when temperature tends to zero, we get 

< i > ~ ((/^o, iy^o) (T^o). 

For the approximate Hamiltonian, one has 

pf ^jfVV _ j^app app j^app _ f app fr app 

n.app'iJk ~ ■'^k ^k ' -^0 ~ yro i-'^app'i^Q 

Similarly to the case of the exact Hamiltonian, 
and 

< A >app^ (^r, ^^r) (T ^ 0) . 

In that way, the Gibbs-Bogolubov inequality at T = reduces to the inequality 



Eo<(^r,^^rj , (110) 

which is called the Rayleigh-Ritz inequality or the variational principle. 

4.8 Bogolubov inequalities 

Bogolubov has proved several more useful inequalities [30lll8lll9]. One of them concerns 
the statistical average of the product of two operators < AB >, with a Hamiltonian H. 
Generally, this average < AB > can be treated as a bilinear form, that is a form linear in 
both A and B. The required conditions are (< A'^A >) > and < Ab >*=< 5+A+ >. 
For statistical averages, these conditions are always valid. Note that A'^A is a semi-positive 
operator. 

Bogolubov inequality for products. The bilinear form < AB >, satisfies the in- 
equality 

<Ab> < < ii+ >< B+B > . (Ill) 

Proof. Let us consider the form 

< (xB^ + yA) (x*B + y*AA > > , 
which is non-negative for arbitrary complex numbers x and y. This can be rewritten as 
\x\^ < B+B > +xy* < 5+i+ > +x*y < Ab > +\y\^ < AA+ > > . 

Setting 

X = - < Ab > , y = < B^B > 

yields 

< AA-^ >< B^B >^ - <Ab > < B^B > > . 
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If < B^B >y^ 0, then the Bogolubov inequahty follows. 

And if < I3^B >= 0, then we set a; = z < AB > and y = 1, with z being a real-valued 
number. This gives 



2z 



<AB> 



+ < AA+ > > . 



Let us send z — )■ — oo. Then the left-hand side tends to — oo, if < AB > is nonzero. This 
is impossible. Hence < AB > = when < B^B > = 0. The Bogolubov inequality again 
holds true. D 

Another inequality has been proved by Bogolubov for the average of the commutators of 
two operators A and B. 

Bogolubov inequality for commutators. For any two operators A and B the in- 
equality holds: 

<[A,B]>'^ <^<A+A + AA-^ > <[B+,[H,B]]> . (112) 

The proof of this inequality is a bit cumbersome, because of which it is omitted here [30] . 

4.9 Nernst law 

Nernst formulated the law as an empirical fact valid for any equilibrium system. But in 
quantum statistical mechanics, this law can be proved. 

Let Co > 1 be the number of quantum states (po possessing the same energy Eq = 
{ipo, Hipo) at T = 0. This parameter (o is called the degeneracy factor. It is independent 
from thermodynamic variables, being defined by the quantum properties of the system. 

Nernst law. The entropy 5* = — Tr p In p for an equilibrium system at zero-temperature 
is 

lim5 = lnCo>0. (113) 

Proof. The statistical operator of an equilibrium system has the form p = e~^^/Tre^'^^. 
When temperature tends to zero, we have 

Tr e"^^ ~ Coe-^^° (T ^ 0) . 

Then limit (113) follows. D 

The quantity 5*0 = In (q has the meaning of the Boltzmann entropy. Since the entropy S 
is an extensive quantity, one has limjv_^oo S/N < oo. Assuming that the degeneracy is finite, 
Co < c>o, the quantity 5*0 can be neglected in the thermodynamic limit when 

lim — = lim = . 

N-^oo N N^oo N 

But 5*0 cannot be neglected if Co oc e"^, with a > 0. 

From the Nernst law it follows that the specific heat at zero temperature is zero. 
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Zero-temperature specific heat. The hmit of the specific heat at zero temperature 

is zero: 

lim Cv = 0. (114) 

r->o 

Proof. At temperature tending to zero, one has 

5*0 = lim S = lim — — = lim — — -— . 

T^O T^O T T^O dT 

Using the expression of the specific heat 

T rdS\ 1 fdE 



we come to the equality 



^"^ N\dT)y N\dT,y 

So = So + N lim Cv . 

From here, limit (114) follows. D 

It is important to emphasize that, for any degeneracy factor C^q, the limits 

lim F = Eq , lim < A > = (ipo, Audq) 

are valid, since Co does not depend on T, and 

lim T In Co = . 

Knowing the behavior of the free energy at low temperature makes it possible to derive an 
equation for the free energy at arbitrary temperatures. To this end, because of the relations 

dF 
F = E-TS, S = - — , 

we have the limits 

OF 
hm F = Eq , hm -— = -^o , 

T-5>0 T-5>0 OT 

serving as boundary conditions for the equation 

dF 
T ^— = F-E . 
dT 

Integrating the latter equation over temperature between zero and T, we obtain the expres- 
sion 

F = Eo-TSo-tJ ^^^dT. (115) 

If 5*0 can be neglected in the thermodynamic limit, then the expression for the free energy 
at large A^ ^ 1 becomes 

F = Eo-tJ ^^^dT. (116) 

The integral over T is well defined because the energy i?, as T — )> 0, tends to Eq not slower 
than as T^. 
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4.10 Parameter variation 

It often happens that the Hamiltonian H = H{X) depends on a parameter A, and one needs 
to find the variation of an observable quantity < A > = Tr pA that is the average of an 
operator A = A{X), which can also depend on this parameter. When the Hamiltonian is a 
function of a parameter, the statistical operator is the function 

P=— - = P(A) 

Tre-/^^ 

of this parameter as well. 

Before considering the variation of the parameter A, let us define the covariance of two 
operators A and B: 

cov(i, B) = - < [A, B]+ > - < A >< B > , 

where the notation for the anticommutator [A, -B]+ = AB + BA is used. The covariance is 
symmetric, 

cov{A, B) = cov{B, A) , 

and additive, 

cov(i + B,C) = cov(i, C) + cov{B, C) . 

The diagonal covariance is the variance, i.e., the dispersion that for a self-adjoint operator 
is 

cov(i, A) = < A^ > - < A>'^ = var(i) . 

For a function of an operator / = f{A), the derivative with respect to the parameter is 
defined as the symmetrized derivative 



M ^(^' ^ i 



9/ dA 



Direct differentiation yields 



Considering the zero-temperature limit, we have the following properties. 
Proposition 1. Let A be a self-adjoint integral of motion, then 

lim /3 cov I A, -- 1 = . (118) 

T^o \ dX 



Proof. A self-adjoint operator is such for which {ipk, Aipp) = {Aip^, ipp) for any ^9^, ipp from 
the Fock space. An operator is the integral of motion when [A, H] = 0. In that case, it enjoys 
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a common set of eigenfunctions with the Hamiltonian, so that Aipk = ^k^k and H(pk = Ek(pk- 
For the normahzed functions, such that {ipk,<^k) = 1) one has d{ipk,(pk)/d\ = 0. Therefore 

d / ;; \ / dA 

— \^(fk, Aipkj = ['^k, -g^ Vk 

When T — )■ 0, then hmx-yo < A> = f(/9o,Ay9o], with y^o being the ground state. In that 
way, 

hm -— < A > = hm < — — > . 

T-i>0 d\ T^o oX 

Then, in view of equahty (117), hmit (118) follows. D 

Proposition 2. For the derivative of the entropy S with respect to a parameter A, the 

zero-temperature limit 

(99 
limT — = (119) 

T^o dX 

is valid. 

Proof. Let us take in Eq. (117) A = H. Then for E =< H >, aX any T > 0, we have 

dE OH ^ f ^y dH 



For the normalized functions. 



d ( fj \ dH 



According to property (118), one has 



lim /3 cov I i/, — - 1 = . 

T^o \ dX ' 



It is straightforward to check that 

dH dE 



<ix>=irx- ^'-T^^^^-'" 



This leads to the relation 



dE dE ^ 1^^ dH 
5A = 5A -^^^M'^' 51 



Since E = E — TS, we get 



d£_d^_rj.dS_ aE _ dH_ dS_ 

'dX~'dX ' dX' 'dX ~'^~dX^^ dX 
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Comparing the above formulas, we obtain 

dS „2 If, SH 

which, according to Eq. (118), gives hmit (119). 

Note that the same result can be derived from the expression S = —Tip In p, from which 

dS ( dH 

Substituting here Inp = /3(F — H) and taking into account Eq. (118) yields limit (119). D 

Proposition 3. Let the Hamiltonian be the sum H = K + W oi the kinetic-energy 
operator and potential-energy operator, respectively. Then the variation of the internal 
energy E with respect to mass m, at zero temperature, gives 

dE 1 ^ dW 

^ = __<i^> + <^> (T = 0). (120) 

om m om 

Proof. Let us set \ = m. From the form of the kinetic-energy operator, we have 

dk 1 ^y dH 1 ^y dW 
K , -— = K + 



This yields 

which is equivalent to 



dm m dm m dm 



dE 1 ^y dW ^ dS 

dm, m dm, dm, 



dE 1 ^y dW 

dm, m dm. 



According to Proposition 2, we have 

ds___ 2 ( fr ^h\ r — -Q 

dm, \ ' dm, J ' t^o dm, 

Thus we come to Eq. (120). D 

Equation (120) shows that the Landau-Lifshitz [28] relation < K > = —m dE/dm, at 
T = 0, is valid only when dW/dm, = 0. However the latter is not always true. For example, 
in the case of a dilute system of cold atoms, one has dW /dm, = —W/m. Then 

-— = H , -— = — var(/f) . 

dm, m, dm, m 

In view of the expression for the specific heat 

1 fdE\ var{H) 



'^^ N \dTjy NT^ 
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we get 


dS 

dm 


^ r 


hm TCv = . 


Therefore, 


instead of the Landau-Lifshit2 


', relat 


ion, we 


! have 






dE 
dm 


E 
m 


(T = 


0). 



(121) 

This shows that one has to be cautious when invoking equations presumed to be well known, 
not forgetting to check the region of applicability of these equations. 

4.11 Kubo formula 

When the system is subject to an external perturbation, its observable quantities change. 
Such a change can be described by the Kubo formula [5"7i[58]. 

Let the system, at time t = 0, be characterized by a Hamiltonian H, and a statistical 
operator p{t). Generally, if is a grand Hamiltonian. In the case when there are no external 
perturbations, the statistical operator satisfies the Liouville equation 

The statistical operator evolution can be represented as p(t) = U{t)p{to)U^{t) by means of 
the evolution operator satisfying the Schrodinger equation 

1 1 m = Hu{t) , 

with the initial condition U{to) = 1. The observable quantity, corresponding to an operator 
A(t) = tJ~^(t)A(to)tJ{t), when there is no perturbation, is given by the average 

< A{t) > = Trp(t)i(to) = Trp(to)i(t) • 

But let us assume that at the moment t > to, the Hamiltonian suddenly changes to Hout- 
Hence now the statistical operator becomes Pout{t) that satisfies the Liouville equation 

i -17 Poutit) = [Hout, Pout{t)] , 

with the initial condition Poutito) = p(to)- Now, the observable quantity, corresponding to 
the operator A, reads as 

< A(t) >out = Ti Poutit) A(to) . 

The change of the observable quantity, called increment, is defined as 

A < A{t) > = < A{t) >out - < A{t) > ■ (122) 

Let the new Hamiltonian have the form Hout = H + AH. Then the Liouville equation 
becomes 

i -Tl Poutit) = [H,pout{t)] + [AH, Poutit)] . 
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Looking for the solution of the Liouville equation in the form 

Poutit) = U{t)D{t)U+{t) 
and introducing the notation 

AH it) = U^{t)AHU{t) , 

we obtain the equation 

Zj^D{t) = [AH{t),b{t)], 

with the initial condition D(to) = Pout{to) = p(to)- The above equation can be rewritten as 
the integral equation 

D{t) = p(to) - I / [AH{ti),Dih)] dti . (123) 

J to 

Using this, the new statistical operator can be represented as 

Poutit) = m-^U{t) I [AH{ti),D{U)] dh U^{t) . (124) 

J to 

Then increment (122) takes the form 

/•* 
A < A{t) >= -i Ti[A{t),AH{ti)]D{ti) dti . (125) 

Jto 

This expression is exact for any AH. 

By constructing an iterative solution for D{t), we start with the zero approximation 
Do{t) = p{to). Then the first approximation is 

Diit) = Doit) -I f [AiJ(ti),p(to)] dti ■ 

Jto 

The second iteration gives 

D2{t) = Di{t) + i-t)^ f r [AH{ti),[AH{t2),p{to)]]dtidt2. 

Jto Jto 

Continuing the iteration procedure, we come to the formal solution 

"tn-l 



D{t)=pito) + J2i-^T dti dt2... 

1 Jto Jto Jto 



n=l J ^0 J^o Jto 



X [AH{ti), [AH{t2), . . . [AH{tn),p{to)]] ■ ■ ■] ■ 
Then increment (125) becomes 

/"* r 1 °° 

A<A{t)>= -i < A{t),AH{ti) > dti + Y^ 

Jto '- -' 1 



t 


rti 


rtu 


)^Hr' / dti 


/ dt2.. 


■ / dtn+lX 


tl Jto 


Jto 


Jto 

(126) 
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X < 



i(t),AiJ(ti) ,AH{t2) ,...,Ai7(t„+i 



> 



If the Hamiltonian perturbation AH = 6H is small, such that 6H — > 0, then one can 
limit himself/herself by the first-order approximation for increment (126). Then A < A >= 
6 < A > is termed the linear response or the linear increment that is 



6<A{t)>= -i < A{t),5H{ti^ 

'to 



> dti 



(127) 



which is the Kuho formula. 

A particular important case is when to — t- — oo, which is called the adiahatic switching. 
In that case, the Kubo formula (127) can be represented as 



5 < A{t) > = -i 
where 6(t) is a unit step function 



e(t) 



0(t-ti) < A{t),5H{t^ 



1, t > 
0, t < 



> dti 



:i28i 



Let us introduce the commutator Green function 



C 



A{t),B{t') =-ieit-t')< A{t),B{t') 



> 



(129) 



which sometimes is named the Bogolubov-Tyablikov Green function. Then the linear re- 
sponse (127) reads as 



5 < A{t) > 



C 



A{t),5H{ti 



dti ■ 



(130) 



If one assumes that the system at t < to has been in equilibrium, with a Hamiltonian H 
and a statistical operator p, then the commutator Green function (129) depends only on the 
difference of times t — t'. This makes it possible to define the Fourier transform 



C 



Ait), Bit') 



X{A,B,uj)e 



'iuj(t-t') 



t') duJ 



27r 



Respectively, the inverse Fourier transform, called the admittance, reads as 



x{AB,u 



C 



i(t),5(0) 



e^^* dt , 



(131) 



where the equality 



e*^* dt = 27r5iuj) 



is used. 

Finally, we come to the expression for the linear response 

6 < Ait) > = lim xiA SH, u) . 
Here all averages are calculated for an equilibrium system. 



(132) 
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4.12 Quasi-equilibrium states 



In section 4.2, it has been mentioned that, using the local information functional, one can 
define a quasi-equilibrium, or locally equilibrium, statistical operator. Here, we shall specify 
the general properties of the corresponding quasi-equilibrium system |46j . 

For any system, the energy Hamiltonian and the number-of-particle operator can be 
written as integrals of the related densities. 



H = / H{x) dx 



N 



I 



h{x) dx 



where h{x) = ■ip^{x)ip{x) is the density-of-particle operator. One can introduce the local 
energy density E{x) =< H{x) > and the local particle density p{x) =< n{x) >. Their 
integrals define the internal energy and the number of particles. 



E = E{x) dx , N = p{x) dx . 
The local information functional can be written in the form 



noc\ 



Trplnp + Ao (Trp — 1) + / Xi{x) TipH{x) — E{x 



dx+ 



+ / ^2(0;) \Tif)h{x) — p{x)\ dx , 

with the Lagrange multipliers guaranteeing the normalization condition and the definitions 
of the local energy density and local particle density. These multipliers can be rewritten 
as Ao = InZ — 1, Ai(s) = /3(x), and A2(a;) = —[3{x)p{x), where /3(x) is the local inverse 
temperature, such that /3(x)T(x) = 1, with T{x) being the local temperature, and p{x) is 
the local chemical potential. 

From the principle of minimal information, we find the local statistical operator 



P 



exp s — / [i{x) H{x) — p{x)n{x) 



dx 



(133) 



with the local partition function 



Trexp 



The quasi-equilibrium potential Q 
potential Q{x), so that 



f3{x) H{x) — p,{x)h{x) 



dx 



In Z can be represented through the local grand 



Q = l3{x)il{x) dx , Z = exp 



f3{x)n{x) dx 



(134) 



From Q = —InZ, we have 



5Q 

S/3{x) 



E[x) — p{x)p(x) 
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While from the first relation with Q{x) in Eq. (134), we get 

— il[x) + p[x) 



Sl3{x) ' ' "^^ ' 6(5{x) 
Comparing the above formulas, we find 

5Vl{x) 



(5/3 (x) 
Because of the equality 



T{x) [E{x) — fi{x)p{x) — Q{x)] . 



6n{x) _ _^, 6n{x) 



5(3{x) ' ' 5T{x) ' 

we obtain 

^^=^ix)[nix)+f^ix)pix)-E{x)] . 

01 [X) 

The local entropy is introduced by the definition 

6n{x) 



Six) 



5T{x) ' 
which results in the expression 

S{x) = /3{x) [E{x) - fi{x)p{x) - n{x)] . (135) 

From here, we have the relations 

^, , 6S(x) , , ^, , 6S(x) 

P{x) = —^ , fi{x) = -T{x) ^ ' 



6E{x) ' 6p{x) 

The Gibbs entropy (information entropy) can be represented as an integral of the local 
entropy: 

S = — Trplnp = / S{x) dx . 
Here the local entropy is as defined in Eq. (135), since 

S = l3{x) [E{x) — p{x)p{x)] dx — Q . 



Taking into account the equalities 



6p{x) ' 6p{x) 6p{x) ' 6p{x) 

we come to the local Gibbs-Duhem relation 

6n{x) = -S{x)6T{x) - p{x)6p{x) . (136) 

Introducing the local pressure P{x) = —Q{x) yields 

6P{x) = S{x)6T{x) + p{x)6p{x) . (137) 

The local quantities naturally appear in nonuniform systems, such as the systems of 
trapped atoms. 
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4.13 Gauge symmetry 

One of the most important properties of systems is gauge symmetry. When the system is 
gauge symmetric, there can be neither Bose-Einstein condensation nor superconductivity. 
The typical energy Hamiltonian is 



H = / ^\x) 



2 



SS + ^W 



iIj{x) dx+ 



+ - / ^lj\x)^lj\x')^{x,x')tp{x')tlj{x) dxdx' , (138) 

where U{x) is an external potential and $(a;, x') is a pair interaction potential. 
The number-of-particle operator is an integral of motion, 

[iV, /J] = , N = I i)\x)i){x) dx . 



This is termed the number-of-particle conservation. 

The global gauge transformation is accomplished by the action of the operator 

U^ = e*'^* , (139) 

where ip is the global phase, being any real number. 

Since [A^, H] = 0, Hamiltonian (138) is gauge symmetric: 

U;HU^ = H . 

In other words, H is invariant under the field-operator transformation 

ij{x) ^ U^i'ixW^ . 

The gauge transformation U^ enjoys the following properties: 

It also possesses the property of additivity: U^pU^t = f/;p+<^'. Therefore the family {U,f} forms 
a group, called the gauge group. 

Among various averages of operator products, one distinguishes two principally differ- 
ent types of the averages. One type, called normal averages consists of the products con- 
taining equal numbers of creation and destruction field operators, e.g., as in the average 
< iIj'^{x)iIj{x') >. The other type of the averages, termed anomalous averages consists of the 
products with unequal numbers of the creation and destruction operators, for instance, as 
in the averages < ipi^x) > or < ■iIj{x)i/j{x') >. 

Statement 1. If the number-of-particle operator is an integral of motion, then all 
anomalous averages are identically zero. 

Proof. The particle-number operator has the following properties of commutation with 
the field operators: 

[ipix), N] = tlj{x) , [ip^x], N] = -^Ij^{x) . 
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From here, one gets 

iIj{x)N = (1 + N)^{x) , ^(x)iV" = (1 + N)"'^{x) . 
The gauge-transformation operator (139) can be represented as 

n=— 

Using this, we come to the relations 

^{x% = e'^^U^^ix) , U^i^{x)lJ^ = e'^ijix) , U^^\x% = e-'^^^^x) . 
In equihbrium, [p, H] = 0. Since N is an integral of motion, it follows that 

[iV, p] = , [f/^, p] = . 
Because of the above equalities, we obtain 

< iIj{x)i/j{x') > = TtpiIj{x)i/j{x') = TtU^pU^iIj{x)iP{x') = 

= 1i pJJ^ il){x)il){x')JJ ^ = e^^'^1ipil){x)il){x') . 

Therefore the equality 

< ^{x)^{x') > (1 - e^'^P) = 

is valid for any if. This implies that 

< ip{x)i;{x') > = , (140) 

while < '?/'^(x)-?/'(x') > 7^ 0. In the same way, it is straightforward to show that 

< tpix) > (1 - e*'^) = 
foe any (p, hence < ipix) >= 0- Generally, we find the equality 

n m 

telling us that it is zero for m ^ n. Because this is zero for an arbitrary ip, we come to the 
conclusion that 

n m 

< JJ?/'^(xi) JJ^(xj) > = {m^n). (141) 

Thus, all anomalous averages are identically zero. D 
Contrary to this, the normal averages are not zero. 



n 



<\{i;\xm<)>^^- 



i=l 



Since the appearance of the anomalous averages necessarily accompanies Bose-Einstein 
condensation or superconducting transition [11112], this means that these transitions cannot 
occur in a gauge-symmetric system, but the symmetry breaking must happen. 
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4.14 Translation invariance 

When the operator of momentum is an integral of motion, there exist constraints on the 
averages, from which it follows that such a system has to be uniform. In this section, we 
understand the variable a; as a spatial variable, x — ?> r. 

Let us assume that there is no external potential, U{r) = 0. Then the operator of 
momentum is an integral of motion: 

[P, H]=0 , P = f ^\r){-iV)^{r) dr . 

This is called the momentum conservation. The family of the translation transformations 

f (r) = exp(iP • r) (142) 

forms the translation group {T(r)}. Recall that, according to Sec. 3.8, one has '^(r) = 
T"'"(r)-?/'(0)T(r). When there is no external potential, the momentum operator commutes 
with Hamiltonian (138), so that 

[f (r), i^] = , [f (r), p] = , 

for an equilibrium p. 

Statement 2. When the operator of momentum is an integral of motion, then 



n m n 



<n^^(r. + r) \{i,{v^+v)> = <Y[i;\v,) \{i;{v^)>, (143) 

for any integers m, n. 

Proof. Starting with the case of two field operators, we have 

< '^/'''"(ri + r)-?/'(r2 + r) > = TrpT^(r)-(/'^(ri)-?/'(r2)T(r) = 

= TrT'^(r)p-?/'^(ri)-?/'(r2)T(r) = Trp-(/'^(ri)-(/'(r2) . 

That is, 

< ■?/'^(ri + r)-?/'(r2 + r) > = < ■?/'^(ri)-?/'(r2) > . 

Similarly, for the product of any number of field operators, we obtain Eq. (143). D 
In particular, we have 

< V'Hri)^(r2) > = < ^^(ri - r2)^(0) > , 

which leads to the uniform particle density 

< ij\r)ij{r) > = < ^\0)ip{0) > . 

For a uniform system of A^ atoms in volume V, it is convenient to use the Fourier basis 
{(^fc(r)} of plane waves 



Vk[r) 



1 Akr 



'V 
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The basis is orthonormalized, {(fk, <^p) = Skp, and complete, 

k 

Imposing the periodic boundary conditions fki'"^ + L) = v^fc(r), where L = {L"} and 
V = Yla -^") on has k ■ L = 27r, hence 

k={r}, A;" = ^n" (n" = 0,±l,±2,...). 

In the thermodynamic hmit, when A^ — )■ oo,^ — )> oo,L° — > oo, the difference between the 
neighbor wave numbers A/c" = 27r/L", for which An" = 1, tends to zero, AA;" — )■ 0. In that 
case, the summation over the wave numbers A;", which is equivalent to the summation over 
n", so that J2k — S{n«}' '^^^ t)e replaced by the integration 



Since 



J2^[dn (n=K}). 



an = L — — , dn = V 



27r ' (27r)3 ' 

the replacement takes the form 



' (2-) 



3 



Passing to the momentum representation, one expands the field operators over the Fourier 
basis, 

^(^) = ^akVk{r) , 

k 

where 

ak = i^k,^) = / (p*k{r)i!{r) dv . 



From the commutation relations ['?/'(r, ipir')]^ = and [ip{r, ■?/'''" (r')]=(i = 6{r — r'), one gets 
the commutation relations 

[ak, Opjzp = , [ak, a],]zp = 6kp • 

Because of the momentum conservation, the system is translation invariant, such that 

d 

— < V^'^(ri + r)V^(r2 + r) > = . 
or 

Passing to the momentum representation in the average 

< ^"^(ri + r)^(r2 + r) > = ^ < ^Op > ipl{ri + r)(pp{v2 + r) 

kp 



yields 
dr 

kp 
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d 

— < ij^ri + r)ij{r2 + r) > = -i ^ < 4ap > (k - p)v?fc(ri + r)v9p(r2 + r) = 

kp 



From here, it is seen that for a uniform system, one has the properties 

< a\ap > = 5kp< a\ak > , < akttp > = 6_kp < a-k^k > ■ (144) 

In general, the average 

n m 

<n< n«p.>^o (145) 

is not zero only when ki + . . . + k„ = pi + . . . + p^- 

In this way, if the operator of momentum is an integral of motion, the system is uniform 
and cannot characterize nonuniform systems, such as crystalline matter. To be able to 
describe the latter, translation symmetry breaking is necessary. 

4.15 Symmetry breaking 

As is stressed in the previous two sections, in order to describe some thermodynamic phases, 
it is necessary to consider systems with broken symmetries. For instance, if the Hamiltonian 
of a magnetic system is invariant under spin inversion or rotation, that is, exhibits the related 
symmetry with respect to the magnetic moment M, then the latter is identically zero [49]. 
Or, when the operator of momentum is an integral of motion, [P, H] = 0, the sys- 
tem enjoys translation symmetry, as is shown in Sec. 4.14, because of which the average 
{ip'^{T)ip{r)) = {ip^{0)ip{0)) does not depend on the spatial variable. But the average 

p(r) =< ^^(r)^(r) > (146) 

represents the spatial distribution of particles, that is, the particle density. If p(r) = p(0) = 
const, then there can be no crystalline matter. 

Similarly, as is discussed in Sec. 4.13, when the system is gauge-symmetric, all anomalous 
averages are identically zero, hence there can be neither Bose-Einstein condensation nor 
superconductivity. 

The above examples demonstrate that, for characterizing some thermodynamic phases, 
one has to look for states with a broken symmetry. The general idea of how the symmetry, 
prescribed by the Hamiltonian, can be broken is straightforward. One has to impose some 
additional conditions, boundary conditions, initial conditions, or like that, for selecting from 
the total Fock space the microscopic states possessing the required symmetry properties. 
This is analogous to the selection of required solutions of nonlinear equations possessing 
multiple solutions. 

The direct way of selecting the required states is by the method of restricted trace. One 
defines a subspace of the Fock space J-'iy G J-", whose states enjoy the necessary symmetry 
properties. And, calculating the observable quantities, one takes the trace only over the 
states of the selected subspace, 

< A >,= Tr^^ pA . (147) 

A general method of symmetry breaking has been developed by Bogolubov |19] and is 
called the method of Bogolubov quasiaverages. The idea of the method is as follows. Let the 
system Hamiltonian H possess some symmetry related to an integral of motion C, so that 
[C, ^] = 0. Suppose we wish to break the corresponding symmetry. For this purpose, we 
introduce the Hamiltonian 

H,{e) = H + eT,, (148) 
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adding to H a, term T^, called source, that breaks the symmetry because it does not commute 
with the integral of motion C, 

[C, f,] ^ , [C, H^is)] ^ . 

The average of the source < F,^ > is assumed to be an extensive quantity proportional to A^ 
for large N ^ 1. The index z/ enumerates the admissible phases. 
The quasiaverage of an operator A{r) is 

< i(r) >,= lim lim < i(r) >^^(,) , (149) 

where 

Trie-'^-^ 
< A>H = ^ • 

Tre-/'^ 

Here limjv_j,oo implies the thermodynamic limit, when A^ — )■ oo, K — )• oo and N jV — ?■ const. 
It is principal that the limits e — > and A^ — t- oo do not commute. 

Instead of the double limit, as in the Bogolubov method, it is possible to employ a single 
limit, as in the method of thermodynamic quasiaverages ^59j, by introducing the Hamiltonian 

H, = H+^t, (0<a<l), (150) 

where the restriction < a < 1 guarantees that the factor l/N"', in the thermodynamic limit, 
goes to zero slower than 1/N, which is necessary for the efficient selection of the states with 
broken symmetry. The average of the source < F^ > is again assumed to be an extensive 
quantity proportional to N for large A^ ::^ 1. Then the thermodynamic quasiaverage is 
defined by the sole thermodjTiamic limit 

< i(r) >^= lim < i(r) >h, . (151) 

For example, to break the spin rotational symmetry, it is sufficient to add to the Hamil- 
tonian the source representing the action of an external magnetic field. To break the trans- 
lation symmetry, which is necessary for describing crystalline matter, it is sufficient to add 
the source 

V,= I ^\v)U{v)^(v)dv, (152) 

in which the external potential is periodic over the lattice vectors {a}, such that f/(r + a) = 
U{v) and [P, F,^] ^ 0. Then the particle density, defined as the Bogolubov quasiaverage 

p^(r) = lim lim < i)\v)i){r) >H,{e) , (153) 

where the index v corresponds to the chosen crystalline lattice vectors {a}, is periodic over 
the lattice vectors, 

p,(r + a) = p,(r). (154) 

Or one can use the thermodynamic quasiaverages (151) defining the single limit 

p^(r) = lim < V"^(r)V'(r) >h, . (155) 
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The method of breaking the translation symmetry by adding a source with an exter- 
nal potential that is periodic over the considered crystalline lattice was, first, advanced by 
Kirkwood [47j. Bogolubov [3^ extended this method to the case of an arbitrary symmetry 
and gave it a precise mathematical meaning. The replacement of the double limit in the 
Bogolubov method by a single thermodynamic limit was suggested in Ref. [59]. Some other 
methods of symmetry breaking are reviewed in Ref. [33] . 

Generally, there can exist several equilibrium states of a system, under the same thermo- 
dynamic variables. However, not all those states are stable. The system always chooses a 
state that is the most stable. When it happens that a state with a broken symmetry is more 
stable than the symmetric state, the system passes to that more stable state with broken 
symmetry. This is termed spontaneous symmetry breaking. 

5 Quadratic Hamiltonians 

Averages can be calculated exactly, when the system Hamiltonian is quadratic with respect 
to field operators. However, such a trivial situation happens rather rarely, corresponding 
to ideal gases. Interacting systems, as a rule, do not allow for exact calculations. But 
approximate calculations are possible by reducing the exact Hamiltonian to an approximate 
quadratic form. 

5.1 Occupation numbers 

In many cases, it is convenient to invoke the representation of occupation numbers. This rep- 
resentation is introduced as follows. Let us, first, pass to the quantum-number representation 
by expanding the field operator over an orthonormal basis, 

i'ix) = ^afcV2fc(a;) , (156) 

k 

where the operators ak and a], generate the Fock space J-'(afc), with k being a set of quantum 
numbers (multi-index). Depending on the particle statistics, the commutation relations are 

[ak, ap]=p = , [ak, al]^ = 6kp , 

the upper sign corresponding to Bose statistics, while the lower, to Fermi statistics. 
The occupation-number operator is 

hk = alttk . (157) 

This is a self-adjoint operator, hi = hk, and it is semipositive, since, for any / G J^{ak), one 
has 

if,nkf) = if,alakf) = {akf,akf) = ||afc/||^ > . 

The number-of-particle operator can be expressed as the sum 

N = ■ip\x)ijj{x) dx = '^ alak = ^ hk 
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of the occupation-number operators. 

For any type of statistics, the commutation relations are vahd: 

[fik, hp] = , [ttk, hp] = Skpttk , [al, hp] = Skpttl , 

where all square brackets represent commutators. 

The eigenvalues z/^ of fik, called the occupation numbers, are given by the eigenproblem 

nk\k >= Vk\k > ■ (158) 

They are non-negative, Uk > 0, since n^ is semipositive. 

For Fermi statistics, the eigenvalues are easily found. As far as afcOfc = —aka^, one has 
al = 0. Then 

nl = alatalak = 4(1 - alak)ak = alat = ht . 

This gives v], = u^, since n1\k >= v1\k >. From here it follows that u^ = 0,1, for Fermi 
statistics. 

But in the case of Bose statistics, the situation is a bit more involved. First, we notice 
that if \k > is an eigenvector of n^, then (afe)™|A; > and (a|.)™'|A; > are also its eigenvectors. 
To prove this, let us consider the equalities 

hkak\k > = {ttkhk + [fik, ttk]) \k> = (z/^ - l)ak\k > , 

hka\\k > = (alfik + [fik, 4] j |/i; > = (z/^ + l)4l^ > • 
Continuing further, we find 

hkiak)"'\k > = (z/fc - m){ak)"'\k > , nkial)'^\k > = (z/^ + m){al)"'\k > . 

The eigenvalues of fik cannot be negative, which means that, after some m, the null vector 
|0 > must exist, such that 

(afc)""|A;> = |0 > , afc|0>=0. 

For the null vector, the eigenvalue of fik is zero. Hence, the eigenvalues Uk of fik are non- 
negative integers z/^ = 0, 1, 2, . . . ranging between and m. The latter can be set to be 
arbitrarily large. 

The occupation-number representation is convenient for the use in the case of quadratic 
Hamiltonians, such as in the case of ideal gases, or when the Hamiltonian can be approxi- 
mately reduced to the quadratic form. Then any function f{fik) of the occupation-number 
operator, acting on the eigenvectors of the latter, yields f{fik)\k >= f{nk)\k >. 

5.2 Noninteracting particles 

The simplest case of a Hamiltonian, quadratic in the field operators, is that characterizing 
the ideal gases of noninteracting particles. Let us consider the Gibbs grand ensemble, with 
the grand Hamiltonian H = H — nN. For noninteracting particles, one has 



H = 4j^{x) 



2 



— + U{x) - fi 
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^{x) dx , (159) 



where U{x) is an external potential. The field operators can be expanded, as in Eq. (156), 
over the orthonormal and complete basis {ipk{x)} of the solutions to the Schrodinger equation 



2 



^ + f^(-)-/^ 



Vk{x) = Ek (Pk{x) . 

As a result, Hamiltonian (159) reduces to the simple form 

H = ^(£fe - iJ,)alak = ^ ujkolak , (160) 

k k 

with the single-particle spectrum Uk = Sk — /i- 

This Hamiltonian is evidently invariant under the local gauge transformations 

ttk -^ ttke''^'^ , 

where ak is a real-valued number. Therefore all anomalous averages are identically zero, 
such as < akttp > = 0, while for the normal averages, we have < al.ap > = 5kp < alak >■ 
The average of the occupation-number operator n^ = alak gives the occupation-number 
distribution 

nk = <hk> . (161) 

To calculate the distribution n^, we need to find the partition function Z = Tie^^^ . To 
this end, we notice that 

Tr exp I -/3 ^ Wfe^afc ) = Tr J]^ exp {-/Sukfik) = H X^ exp(-/3a;fcZ/fc) , 

\ k / k k Uk 

with z/fc = 0, 1,2, . . . for Bose statistics and z/^ = 0, 1 for Fermi statistics. Hence, for Bose 

statistics, 

E^-fiui^Uk = 1-1- Q-I^^k _|_ g-2/3tjfe I _ 

1 — p^ 



Consequently, 



For Fermi statistics. 
Therefore, 



Z = Y[{l-e-^^'y^ (Bose). 

k 

\~^ Q-I^'^kl'k =1-1- Q-P'^k 



^k 



Z = JJ (l + 6"'^'^'=) (Fermi) . 

k 

Combining both these cases, we get 

Z = JJ(lTe-'^"'=)^' . (162) 

k 

The grand potential is 

fi = -T In Z = -T In Tre"^^ = ±T ^ In (l ^ 6"^'^'=) . 

k 
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From here, we find the occupation-number distribution 

J\0 1 

nk = <hk> = - — = -3- . (163) 

6uJk eP'^k ip 1 

Using the identity e^^'^ = (1 ± nk)/nk, we get the grand potential 

f] = TT^ln(l±nfc) . (164) 

k 

Other thermodynamic quantities can also be calculated exactly using the derivatives 

—— = (3'^ujknk{l±nk) , —— = I3nk{l ± Uk) . 

Thus the internal energy is 

E = < H > +n = ^ EkUk . 

k 



For the entropy, we have 



y^^ k 

The number of particles is 



S" = - ( — j = J^ [PuJkUk ± ln(l ± Hk)] . 






k 

Since fi = —PV, the pressure becomes 

k 

The free energy is F = E — TS = Q + ^N. 

It is worth recalling that fik is a semipositive operator, hence the distribution nk =< fik > 
has also to be semipositive. This implies that e'^'^'= > ±1. For Fermi statistics, n^ > in 
any case. But for Bose statistics, it should be e'^'^* > 1, which means that Uk > 0, so that 

H < Ek. 

The specific heat reads as 

Cy=N{-df] =|rE^^Ml±n.), (165) 

where the derivative 

— = l3^^ujlnk(l±nk) 

k 

is used. For a stable system, it should be: < Cy < oo. 
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For the isothermal compressibihty, using the derivative 

o AT" 

-Q- = P^nk{l ±nk) , 

^ k 



we have 



1 (dN\ var(iV) 



with the number-operator dispersion 

var(iV) = ^nfc(l±n,). (167) 

it 

The system is stable provided that < /ty < oo, which implies that var(A^) ~ A^. 

Similar expressions arise for interacting systems treated in a mean-field approximation, 
with the difference that ook becomes an effective spectrum defined by the corresponding 
self-consistency conditions. 

5.3 Correlation functions 

In the case of the quadratic grand Hamiltonian (160), it is straightforward to calculate 
various correlation functions. Let hi be either ah or al . And let us consider the correlation 
functions of the type < 6162 •• • ^2n >, containing even numbers of the field operators. Such 
correlators can be expressed through the combinations of the lower-order correlators, as we 
shall illustrate for the fourth-order correlation function. 

Wick theorem. The correlation function < 6162&3&4 >j where the averaging is defined 
using the equilibrium statistical operator p = e~^^ /Z, with Hamiltonian (160), reads as 

< 61626364 > = < 6162 >< 6364 > ± < 6163 >< 6264 > + < 6164 >< 6263 > . (168) 

Proof. With Hamiltonian (160), we have 

QkH = {H + Uk)ak , afc-ff" = {H + UkYau , 



from where 
Therefore 



a,e-^^ = e-^^a^e-^-^ , a^e"^^ = e'^^ale'^-^ 



(1 ± nk)akf) = rikpak , Ukolp = p4(l ± Uk) . 
Using this, for any operator A, we get 

(1 ± Uk) < Aak > = Uk < akA > , Uk < Aal > = (1 ± n^) < alA > 

< flfei > = (1 ± Hk) < [flfc, A]^> , < alA > = ^nk< [al, A]^ > 
Recall the commutation relations 

[ofc, Opjzp = , [aj., aj]=p = , 
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and the properties 

< alap > = dkpUk , < akal > = 6kp{l ± n^) . 
Excluding here Skp, we find 

(1 ± nk)[ak, bi]^ = < akk > , =Fr^fc[4, k]^ = T < alk > . 
Using the equahty 

[A, BCD]^ = [A, B]^CD ± B[A, C]^D + BC[A, D]^ , 
vahd for any operators A,B,C,D, we come to the expression 

< akb2hbi > = < akb2 >< 63^4 > ± < ctfc&s >< ^2^4 > + < akb^ >< 6263 > , 

and to the same expression for a|,. This proves the theorem. D 

Because of relation (156), the equivalent representation is valid for the field operators 
'il>{x). Let ipi be either ip{xi) or ^^(xj). Then 

< V'1^2^3^4 > = < ^1^2 >< V^3^4 > ± < IpllpS >< i^2'^4 > + < Ipllpi >< V'2V'3 > • 

Generally, any correlator < 6i&2---^2n >, calculated with the quadratic Hamiltonian 
(160), results in the combination 

< 616263^4 ...b2n> = < &I&2 >< &3^4 > . . . < &2n-l&2n > ± 
± < 6163 >< 62^4 > . . . < &2n-l&2n > + < &1&4 >< &2&3 > • • • < &2n-l&2n > + • ■ ■ ■ 

In the same way, the correlator < ilJiip2 ■ ■ ■ V'2n > can be represented as the combination of 
the products of all possible pair contractions < ipi'ipj >■ 

However, one has to be cautious when the quadratic Hamiltonian is not exact, but is a 
result of an approximation. Then calculating the correlators of the orders higher than second 
can lead to unreasonable conclusions. This is the simple consequence of the approximation 
theory prescribing that, if the main Hamiltonian is of second order, then only the correla- 
tors of second order are reliable. In brief, accepting an approximation of some order, one 
should not go outside of the region of applicability of that approximation, but must restrain 
calculations to the same order. 

5.4 Canonical transformations 

The general form of a quadratic Hamiltonian is 



H 



iIj''{x)A{x,x')iIj{x') + - B{x,x')iIj\x)iIj\x') + - B*{x,x')iIj{x)iIj{x') 

Zi Zi 



xAjtJL/ \JbJu • 



(169) 
Such a quadratic form can be diagonalized by means of canonical transformations [30] . 
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The Hamiltonian can be rewritten in the quantum-number representation, for which 
one expands the field operators over an orthonormal complete basis, labeled by quantum 
multi-indices k, 

k 

Introducing the notation 

A{x, x') = ^ Atp'fkix)ip*{x') , B{x, x') = ^ Bkp^k{x)Vp{x') , 

kp kp 

one gets the Hamiltonian 

H = ^i Akpalttp + - Bkp alal + - Blpttpttk j ■ (170) 

kp ^ ' 

In general, the Hamiltonian can contain nonoperator terms, which are not included here. 

The Hamiltonian is self-adjoint, H^ = H, which requires that Af = Apj.. Also, it should 
be that Bp^ = ^B^p, since 



Yl ^'^P^'Wp = ± 1] Bkpttlal = ±Y1 ^pfc«I« 



t 
p 

kp kp kp 



The general form of the Bogolubov canonical transformations is 

ak = Yl i^kpbp + vlpbl) , al = Y^ [ulj^l + v^php) . (17i; 



p 



It is termed canonical because it does not change the particle statistics, so that the commu- 
tation relations for the operators au are the same as these for the operators hk'. 

[ofc, ap]zp = , [cfe, a],]zp = 5kp , 

[hk, hp]^ = , [bk, 6],]=F = Skp ■ 
This imposes the restrictions 

Yl {ulqVpq T VkqU*pq) = , ^ {ulgUpg T VkqV*pg) = 6kp ■ 

One requires that the Hamiltonian would acquire the diagonal form called the Bogolubov 
Hamiltonian 

H = J2^kblbk + Eo. (172) 

k 

The diagonalization condition gives 



Eo = Eq = — y ^ujk\vpk\ 

kp 
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while the spectrum Uk = ujI and the expansion functions are defined by the Bogolubov 
equations 

y^ [{Akg - 5kgUJp)Uqp + BkqVqp] = , ^ [{^Iq + 5kqUJp)Vqp + Bl^Uqp] = . (173) 

q q 

The operators a^ describe particles, while the operators bk correspond to quasiparticles, 
that is, to excitations. The distribution of particles is 

Uk =< alttk >= ^ [(l^fepl^ ± \vkp\^) TTfc + \vkp\^] , 
p 

while the distribution of excitations is 

TTfc = < blbk > 



(,l3^k ip 1 • 



The particle distribution is summed to the total number of particles A^ = Ylk''^k- But 
the sum over the excitation distribution is not fixed, being dependent on thermodynamic 
variables. 

Because of the diagonal form of the Bogolubov Hamiltonian (172), one has 

< fefe > = , < bkbp > = , < blbp > = SkpTik , < bkbl > = 4p(l ± tt^) . 

As follows from transformation (171), < aj, >= 0, that is, there is no generation of single 
particles. However, there can exist the anomalous average 

< akttp > = Y1 [{vlqUpq ± UkqVpq) TTq + UkqV*q\ . 

q 

The nonzero < akap >y^ implies the breaking of gauge symmetry for pairs, which can be 
created or annihilated together. Respectively, the anomalous average 

< ^{x)^{x') > = 5Z < dkttp > ipk{x)'^pix') 

kp 

is also nonzero. 

5.5 Linear terms 

It may happen that the Hamiltonian contains the terms linear in the field operators, as in 
the expression 

H = J2[ ^kpalttp + - Bkpalal + - Bl^apttk j + J2 {^^al + C^Cfc) . (174) 

kp ^ '' k 

This Hamiltonian is not gauge invariant, because of which < a^ > is, generally, nonzero, 
meaning that these particles can be created and annihilated. Since the particle spin is 
conserved, such a situation can occur only for particles of spin 0, that is, for Bose particles. 
Therefore, the operators a^ here satisfy the Bose commutation relations. 
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To reduce Hainiltonian (174) to a diagonal form, one, first, diagonalizes the quadratic 
part of the Hamiltonian as earher in the previous section, which yields 

H = J2 {''kblbk + Dkbl + D*6fc) + Eo , 

k 

with 

Dk = Y^ {CpU*pk + C*Vpk) . 
p 
Then one invokes another canonical transformation 

bk = bk + — , (175) 

^k 

where the operators b^ satisfy the same Bose commutation relations as the operators a^, 

[bk, bp]=0 , [bk, bl] = 6kp . 
As a result, we come to the diagonal form 

H = J2^kblh + Eo, (176) 

k 

in which 

I n |2 

Eq = Eq — y . 

k ^^ 

The diagonal form of Hamiltonian (176) tells us that 

< 6fc > = , < bkbp > = , 
and the distribution of excitations is given by the equation 

< blk > - ^'' 



^k"P 



e/3'^fc ip 1 ■ 



However the quasiparticles, characterized by the operators 6^, are not conserved and can 
be created or destructed as single objects as well as pairs: 



< 6fc > = - — , <bkbp> = ^^ , < hlh > = < blh > + ^^ 



— , < hbp > = 



The same concerns the original particles described by the operators a^, which can be created 
or annihilated by pairs as well as by single particles. For example. 



Y- fDp d; ^ \ 

< Ofe >= - > Ukp + — %p 

^ \UJp Up J 



This shows that, if the Hamiltonian contains the terms linear in the field operators, as in 
Eq. (174), hence the gauge symmetry is broken, then < a^ > is not zero. 
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5.6 Hartree approximation 

The Hainiltonian of interacting particles is not quadratic and, generally, is given by the 
expression 



H 



■?/'^(x) 






tlj{x) dx -\ — / ilj\x)^lj\x')^{x, x')tlj{x')tlj{x) dxdx . 



(177) 

To reduce it to a quadratic form, one needs to resort to some approximations. The simplest 
such an approximation is the Hartree approximation that, for averages, corresponds to the 
decoupling 

< i)\xi)i)\x2)i){x^)i){xi) > = < i;^{xi)ip{x4) >< i)^ {x2)ij{xs) > . (178) 

In terms of the field operators, this is equivalent to the equality 

ij\xi)^^{x2)^{x3)^{x4) = ^^(xi)V'(a:4) < ^''{X2)^{X3) > + 

+ < ip\xi)ip{x4) > ^\x2)ij{x3)- < ^\xi)ij{x4) >< ij\x2)^p{x3) > . (179) 

Introducing the Hartree potential 

and the average density 

p{x) = < h{x) > , h{x) = tlj\x)^lj{x) 
makes it possible to reduce the initial Hamiltonian (177) to the Hartree approximation 



H 



ip^x) 



-7^ + U{x) + Uh{,x) - n 
2m 



ip{x) dx + Eq , 



(180) 



where 



Eq = / ^{x,x')p{x)p{x') dxdx' 



If we choose the orthonormal basis {ipk{x)} as the set of solutions to the equation 



- -— + U{x) + Uh{x) 
2m 



Vk{x) = Sk^fkix) , 

then, expanding the field operators iIj{x) over this basis, we come to the Hamiltonian 

iJ = ^ Ukalak + Eq {uk = Sk - p) ■ 

k 

For the average density, we have 

P(^) = X]^fc|v'fe(a;)|^ , nfe = < 4afe > = (e^'^'= =F l) 

k 



'i8r 
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As an example, let us consider a uniform system, when there is no external potential, 
U{r) = 0, and let x = r. Also, we take into account that the interaction potential usually 
depends on the difference of the spatial variables, <l>(r, r') = <l>(r — r'). Then the expansion 
basis is composed of the plane waves (fkii^) = e^^'^/y/V and the average density is 

P(r) = y ^nk=p, 

k 

which is a constant defining the chemical potential fi = fi{T,p). Consequently, 

Unir) = p$o , %= $(r) dr , Sk = — + p% ■ 
The internal energy of the system is 

E = <H>+pN = Eo + J2^kek, Eo = - -Np%, 

k 

with Eq being the ground-state energy. 

As is seen, for a uniform system, the Hartree approximation is applicable only if the 
interaction potential is integrable: 



<l>fr) dr 



< oo . 



Another example is a periodic system with a crystalline lattice defined by the crystalline 
vectors {aj}, with i = 1,2,. . . ,N enumerating the lattice sites. Then the expansion basis 
can be formed by well localized functions v^i(r) = (/^(r — a^), called localized orbitals. Good 
localization means that the overlap integral 

(v?i, ^Pj) = / V9*(r - ai)ip{r - a^-) dr 

is small. By the Cauchy-Schwarz inequality, 

\{fi,fj)\'^ < \{Vi,Vi){Vj,Vj)\ ■ 

Since the orbitals are normalized, {ipi,ipi) = 1, one has \{ipi,ipj)\ < 1. And the orbitals are 
well localized, when \{ipi, ipj)\ ^ 1 (i 7^ j). If the orbitals are so well localized that they do 
not intersect in space, then the requirement for the integrability of the interaction potential 
is not compulsory. 

5.7 Hartree-Fock approximation 

The Hartree approximation can work well when the particle exchange is negligible, for in- 
stance, as in the case of well localized orbitals. When the particle exchange is important, 
one has to resort to the Hartree-Fock approximation. This approximation implies for the 
averages the decoupling 
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± < ^^(Xi)^(x3) >< i'^{x2)^{x4) > . (182) 

And in terms of the operators, this is equivalent to the construction 

ij^Xi)^^{x2)^{x3)^{x4) = ^■^(Xi)V'(x4) < ip''{x2)^{x3) > + 

+ < Ip^ {xi)iIj{x4) > llj\x2)lp{x3)- < Ip^ {xi)iIj{x4) >< Ip^ {X2)lp{x3) > ± 
±'llj''{Xi)tlj{x3) < ^P\X2)^P{X4,) > ± < ^\xi)lp{xs) > 1p\x2)ij{x4)T 

=F < ij\xi)i){x3) >< i'\x2)i'ixi) > . (183) 

Let us define the Hartree potential Uh{x), as in the previous section. And let us introduce 
the Hartree-Fock operator potential 

Uhf{x) = Unix) + Uex{x) , 
where Uex{x) is the exchange operator potential defined by the action 

Uexix)^lj{x) = ± / ^{x, x')pi{x, x')tlj{x') dx' . 



Then the Hartree-Fock approximation for Hamiltonian (177) leads to 



H 



^'^(x) 



h U{x) + Uhf{x) - fi 

2m 



ip{x) dx + -E, 



HF 
; 



:i84) 



where 



S, 



HF 



- / $(a;, x') \^p{x)p{x') ± |pi(a;,a;')|^] dxdx' 



and the first-order density matrix is 

pi{x,x') = < ip'^{x')^p{x) >= '^nk'^kix)ipl{x') . 

k 

One can choose the expansion basis made of the functions that are the solutions to the 
equation 

- -— + U{x) + Uhf{x) 



2m 



Vk{x) = Ek^kix) . 



The difference with the Hartree approximation is that now Uhf{x) is a nonlocal potential. 
For a uniform system, when x — )■ r and U{r) = 0, the expansion basis is composed of the 
plane waves. The exchange operator potential acts on the latter according to the rule 

Uexir)M^) = ± /" $(r - r')pi(r, r')v^fe(r') dr' . 

Assuming that the interaction potential is integrable, one invokes the Fourier transform 

$(r) = ^ E -^fce^^-^ , ^k= f $(r)e"^^-^ dr . 
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Hamiltonian (184) becomes 

H = $^(^fe - /^)4«fe + E^^ , (185) 



with the spectrum 

and the nonoperator term 



2m • "--'V /^--^-P 



\ kp 



Eq^ = - - I ivpcpo ± — > ^ij^^^nunp 



2 
In the case of an isotropic system, one has $_fc = $fc and n^k = nk- 

5.8 Hartree-Fock-Bogolubov approximation 

The Hartree-Fock approximation preserves the gauge symmetry. Therefore it cannot be 
apphed for systems where gauge symmetry becomes broken. That is, the Hartree-Fock 
approximation is not apphcable for systems with Bose-Einstein condensate as well as for 
superconducting systems. 

The approximation that takes into account possible gauge symmetry breaking is the 
Hartree-Fock-Bogolubov approximation that for the averages reads as the decoupling 

< ^P\xi)tlj\x2)^p{x3)tlj{x4) > = < ^P\xi)^Ij{x4) >< ^p\x2)^{xs) > ± 
± < ^lj\xi)ilj{x3) >< tp\x2)^{x4) > + < ^''{Xi)^\x2) >< ■?/'(x3)-?/'(x4) > . (186) 

This differs from the Hartree-Fock approximation by the term containing the anomalous 
averages < ■ip{x)i/j{x') > and < ip'^ {x)ip^ {x') >. To express the Hartree-Fock-Bogolubov 
approximation in the operator terms, we shall use the short-hand notation for the field 
operators ^i = ^{xi) and ipj = ^"^{xi). Then this approximation corresponds to the operator 
expression 

i^llpllpSlpi = ^1^4 < ^1^3 > + < ^l^i > ^I'ipS- 
- < iIj\iIJ4 >< V'1^3 > i^lVa < ^1^4 > ± < ^1^3 > ^1^4T 
T < ^1^3 >< ^"1^4 > +^iVI < ^3V^4 > + < ^1^1 > ^"3^4- < ^iVI >< ^3^4 > • (187) 

Hamiltonian (177) reduces to the form 



H= ^\x) 



2 



V 

h U{x) + Uhf{x) - /i 

Im 



ip{x) dx+ 



+ - / ^{x,x') [il)\x)^\x') < i){x')^{x) > +tlj{x')^lj{x) < ^/j\x)tlj\x') >] dxdx' + E^^^ , 

(188) 
in which 

Eq^^s = E^^ - ]- /"$(x,x')| < ij{x')i){x) > |2 dxdx' . 



Expanding the field operators over a basis and introducing the notation 



^kp = / vlix) 



+ U{x) + Uhf{.x) - yu 



(Pp{x) dx , 



Bkp = / ipl{x)^{x,x') < ip{x')ip{x) > (Pp{x') dxdx' 
we can transform Hamiltonian (188) to the expression 

H = ^i Akpalttp + - Bkpttlal + - Bl^apttk \ + E^ 



HFB 




(189) 



kp 



The operator part here is equivalent to that in Eq. (170), hence can be diagonalized in the 
same way as in Sec. 5.4. 

As in the previous cases, the interaction potential, generally, has to be integrable for the 
applicability of this approximation. 

5.9 Kirkwood approximation 

Quite often, the interaction potentials are not integrable, because of which the direct appli- 
cation of the mean-field approximations of the previous sections can be inadmissible. 

An example of a nonintegrable potential is the singular power-law potential $(r) oc 1/r", 
where r = |r|. In the three-dimensional space, the integral 

/"^ dr 
$(r) dr oc 







n-2 



r 



where R is the system radius, diverges at zero, if n > 3. This is classified as ultraviolet 
divergence. 

A rather popular interaction potential is the Lennard- Jones potential 



$(r) = 4e 



0"\12 /(J 



describing the interactions between many neutral particles, such as atoms and molecules. 
The nonintegrable potentials are often termed the hard-core potentials. In the case of the 
Lennard- Jones potential, a is the hard-core radius. 

Kirkwood |17] suggested a method of dealing with hard-core interaction potentials by 
employing the modified Hartree decoupling 

< i)\v)i)\v')i){r')i){v) > = g{r, r) < ^^(r)V'(r) >< ^^(r')^(r') > , (190) 

which differs from the Hartree approximation by the factor g{r, r') that takes into account 
short-range particle correlations, tending to zero, where the interaction potential tends to 
infinity. The actual form of g{r,r') should be defined from other physical assumptions. In 
the operator representation, the Kirkwood approximation reads as 

ij\r)i;\r')ilj{r')4j{r) = g{r,r') [^^(r)^(r) < ^^(r^r') > + 
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+ < 7/'^(r)^(r) > ^^(r')^(r')- < ^^{r)ip{r) >< ^^(r')^(r') >] . (191) 

Substituting form (191) into the Hamiltonian, one meets the effective Kirkwood potential 



UM 



j '^(r,r')5'(r,r')p(r') dr' , 



which is similar to the Hartree potential, but with the principal difference that now, instead 
of the bare interaction potential, one has the smoothed potential 

$(r,r') = ^(r,r')$(r,r') , 

The latter is integrable. 



<l>('r,r') dr' 



< oo 



even if the bare interaction potential is not. 

Equation (190), as is evident, is nothing but the definition of the pair correlation function 
g{r, r'). Kirkwood suggested to extract the pair correlation function g{r, r') from experiment. 
This could be done by measuring the structure factor 



S{k) = ^ f [< n(r)n(r') > -p(r)p(r')] e^^'^^'-'"') drdv' 



(192) 



where n(r) = ■?/'''"(r)-?/'(r) is the density operator. The density- density correlation function 
can be expressed through the pair correlation function: 

< n(r)n(r') > = p{r)6{r - r') + p{r)p{r')g{r, r') , (193) 

which gives 

S{k) = 1 + ^1 p{r)p{r') [g{r, r') - l]e-'^- ('-'■') drdr' . (194) 

In the uniform case, when p(r) = p and g{r, r') = g{r — r'), one has 

S{k) = l + p f[g{r) - l]e-'^-'' dr . (195) 

Taking the inverse Fourier transform yields 

^(r) = l + ^ 5^ [5(k) - 1] e^*^- . 

k 

For a large system, one passes from summation to integration, obtaining 

^(r) = l + iy"[5(k)-l]e^'^-^. (196) 

Measuring in experiment the structure factor S(k), one can calculate the pair correlation 
function by means of relation (196). 



80 



5.10 Correlated approximations 

The Kirkwood idea can be generalized as follows [BUj- For simplicity, we again use the 
notation of Sec. 5.8 denoting the field operators as ipi = i'i'^i)- And let us employ a brief 
notation gu = 5'(ri, r2) for the pair correlation function. The general form of a correlated 
approximation is 



V^i^2^3V^4 = y/9i29u -^m^sip^ , (197) 

L J MF 

where the term in the square brackets in the right-hand side implies some of the mean- 
field approximations. For instance, using for [. . .]mf the Hartree approximation leads to 
the Kirkwood approximation (190). Employing for [. . .]mf the Hartree- Fock approximation 
results in the correlated Hartree-Fock form 

+ < ^^(r)V(r) > V^(r')V(r')- < ^^(r)^(r) >< V^(r')V(r') > ± 

±V^t(r)^(r') < ^t(r')^(r) > ± < ^lj\r)ip{r') > ^^(r')^(r)^ 

T < ^^(r)^(r') >< ^^(r')^(r) >] . (198) 

In the same way, one can introduce the correlated Hartree-Fock-Bogolubov approximation. 
One can remember that the pair interaction potentials depend on the difference between 
the spatial variables, so that $(r, r') = $(r — r'). Similarly, the pair correlation function 
also usually depends on this difference, g{r,r') = g{r — r'). Thus the smoothed interaction 
potential becomes $(r) = g{r)^{r). Since g{r) -^ when $(r) -^ oo, the smoothed potential 
is integrable. 



$(r) dr 



< oo . 



The pair correlation function is real-valued and symmetric, such that g{—r) = g{r). It 
can be expressed through the second-order density matrix as in Sec. 3.10. Thus, for a 
uniform system, it is 

Mr.r-r.r') 

P 

The behavior of the pair correlation function is such that g{r) — )■ 1, when $(r) — > 0. In view 
of the normalization 

f P2ir,r',r,r')drdT' = NiN-l) , 

the property 

^ j g{v)dv = l- 1 (200) 

follows. 

Instead of finding an approximation for the second-order density matrix, thus obtaining 
the pair correlation function (199), one often directly invokes some approximations for the 
latter. The simplest such a way is to use the cutoff approximation 

, . / 0, r < a 
^^ ' 1 1, r > a . 



Sometimes, one resorts to the classical approximation 

g{r) = exp{-/3<l>(r)} , 

which can be used only for high temperatures. 

For quantum systems, it is possible to use the quantum approximation, by setting g{r) = 
\(p{r)\^, with the function ip{r) being defined by the Schrodinger equation for the relative 

motion: 

' V^ 1 

+ $(r) + E 



2mQ 



yp(r) = 



where 

mim2 m 

rriQ = , mi = m2 = m , rrao = — 

nil + ^^2 2 

If at short distance r — )■ 0, the interaction potential behaves as 



a 



r 



12 



$(r)~4£ - (r^O), 



then the solution to the above equation yields 



(^(r) ~ exp <j - ^ (^) \ (r ^ 0) 



with the correlation parameter k = 4vm£?(T^/5. Therefore 

g{r) ^ exp <^ -k (^- 

One also uses variational methods for constructing the pair correlation functions and 
there exist other more complicated ways of finding their approximate expressions [16] . When 
correlations and interactions between particles are strong, it may happen that even integrable 
potentials require to take account of the pair correlation function. 

Different regimes of describing the system can be classified by considering characteristic 
lengths. The interaction radius is given by the expression 

[^ |r|$(r)(ir 
ro = a + T L\^ ■ (201) 



The scattering length 



^ I $(r) dr (202) 



characterizes the intensity of interactions. These lengths should be compared with the mean 
interparticle distance a and the thermal wavelength Xt, 



Particle interactions are strong, when a^ ^ a, while interaction correlations are strong, 
if To ^ a. For At -C a, one has the classical regime. And the quantum regime develops when 
At ~ a, becoming essentially quantum for At > a. 
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The approximations of the mean-field type provide sufficiently accurate description, when 
the particles in the system are strongly correlated, so that the inffuence of each of them 
is extended onto many other particles. In such a case, each particle can be imagined to 
be subject to a mean field created by other particles. Generally, there exist two kinds of 
correlations, interaction correlations and coherence correlations. 

Interaction correlations are strong, when the interaction radius Tq is much larger than the 
interatomic distance a. Then each particle experiences the action of many other particles of 
the system. Such interactions are termed long-range. 

Coherence correlations characterize the particle correlations that develop not merely be- 
cause of particle interactions but also due to particle quantum statistics. The coherence 

radius is defined as 

^ fr\Cir,0)\dr 

^"'- f\C{T,0)\dT ' ^^^^^ 

where the correlation function 

C(r,r') = ^iML (204) 

VP(r)p(r') 

is expressed through the density matrix pi(r, r'). The coherence correlations are strong, if 
the coherence radius Vcoh is niuch larger than the interatomic distance a. 

Summarizing, particle correlations in the system are strong, provided that at least one 
of the conditions holds true, either the interaction correlations or the coherence correlations 
are strong, so that either ro 3> a or Vcoh 3> a. That is, particle correlations are strong, when 
at least one of the following parameters is large: 

^ » 1 , ^ » 1 . (205) 

a a 

In that case, it is admissible to resort to a mean-field type approximation. 

In conclusion, this Tutorial presents the basic notions of quantum statistical physics that 
are necessary for the correct description of quantum atomic systems. Accurately applying 
these notions will help to avoid mistakes that, unfortunately, are rather common in the 
current literature. The exposition in the Tutorial has been aimed to be, from one side, 
sufficiently detailed to be easily read and, from another side, quite brief, omitting lengthy 
discussions that could be found in the cited references. In the following parts, the material 
of this Tutorial will be applied for a more concrete description of cold atoms satisfying 
Bose-Einstein and Fermi-Dirac statistics. 
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